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FOREWORD 


The program reported herein is a continuation of prior 
work under Contract NAS3-20615 from the NASA Lewis 
Research Center reported in NASA CR-159^65 in December 
1978. Both programs were carried out by the same 
principal personnel. However, Barry M. Cohen, who 
was a member of the Comstock & Wescott staff during 
the prior program, has acted as an independent consult- 
ant to Comstock & Wescott during the present program. 
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SUMMARY 


The objective of this work was to extend the prior work 
reported in NASA CR-159465, December 1978, on the de- 
velopment of a phase-change thermal storage system (TESS) 
using modified anhydrous sodium hydroxide and a single, 
passive, tube-intensive heat exchanger, for charging 
and discharging heat at approximately 584 K by a non- 
phase-change heat transfer liquid, Therminol-66 . 

The three parts of the present program were (1) chemical 
analyses on the storage medium used in the experimental 
model of the prior work, and examination of its heat 
exchanger, (2) production of additional experimental 
data using a different non-phase change heat transfer 
fluid ( Caloria HT-43) to extend the validation of a 
computer model of the TESS, and (3) development of a new 
computer model of a TESS containing an additional heat 
exchanger for the vaporization and superheating of a 
power fluid for a Rankine cycle power generator, and the 
use of the model to develop a cost optimized reference 
TESS design for comparison with the single heat ex- 
changer TESS design. 

From the results, it was concluded that during the prior 
experimental work, no chemical degradation of the medium 
occurred, but some physical segregation of components 
was observed which had no apparent deleterious effect on 
the performance of the system. Physical segregation can 
be reversed by fully remelting the medium. No signifi- 
cant deterioration of the heat exchanger was found. 

A computer model was successfully developed which provides 
for 3-way heat transfer, between each fluid and the med- 
ium, and directly between the two fluids. 

The reference design of the 2-heat exchanger TESS (as 
compared to the single heat exchanger TESS) improved 
thermal performance, a lower capital cost (ss. 58 #) when 
credited for the cost of the external heat exchanger 
which it displaces, and improved overall system opera- 
ting cost. 




INTRODUCTION 


Background 


The work reported herein is an extension of prior work 
under two completed contracts: 

1. The development of a computer model (CM) 
of a phase-change thermal energy storage 
system (TESS) under Purchase Order No. 

87-5030 from Sandia Laboratories, Liver- 
more, California, completed in December 
1976. 

2. The development of a phase-change thermal 
storage system using modified anhydrous 
sodium hydroxide for solar electric power 
generation, under Contract NAS3-20615 
from NASA Lewis Research Center, which was 
reported in NASA CR-159^65 in December 
1978. 

This prior work included the design, construction, and 
testing of an experimental scale-model of a TESS suit- 
able for use in mid-temperature range (approximately 
580 K) solar powered electricity generating systems. The 
model was specifically designed for reliable scaling up 
to meet the requirements of the Solar Total Energy Test 
Facility (STETF) operated by the Sandia Laboratories at 
Albuquerque, New Mexico. 

This was used to generate experimental data with which 
to verify the predictions of the CM, which was modified 
in several respects until acceptable agreement was ob- 
tained . 

The TESS studied in the prior work consisted of a tank 
containing the phase-change medium, and a single "tube- 
intensive" heat exchanger, which is used both for charg- 
ing and discharging heat by means of a non-phase-change 
heat transfer liquid, Therminol-66 (T-66). In the 
STETF, heat withdrawn from storage is delivered via 
the T-66 to an external heat exchanger which produces 
superheated toluene vapor to drive a turbine. 
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The nominal composition of the phase-change thermal 
storage medium (Thermkeep*) is: 

Anhydrous NaNO, Commercial Grade 91*8$ Cwt 1 

NaNC> 3 8.0 

Mn0 2 0.2 

The commercial grade of NaOH typically contains 1-2$ of 
NaCl and 1/2-1$ of Na 2 C0 3< 

Phase-change TESS using Thermkeep have been under devel- 
opment by Comstock & Wescott for many years in applica- 
tions using electric heaters and operating at maximum 
temperatures as high as 725 K (see, for example. Ref. 1) . 
However, the referenced work was outside the range of 
previous experience in three respects. First, the heat 
exchanger involved more and smaller tubes which required 
a different method of support. Second, the cycling mode 
of the TESS left a portion of the Thermkeep in the bot- 
tom of the tank unmelted, and because Thermkeep is a non- 
eutectic mixture, this might lead to segregation of its 
components. Third, as explained in detail later, for 
chemical stability in contact with a steel vessel and 
heat exchanger, the system must have contact with air 
to prevent reduction of NaN0 3 to NaN0 2 . This was pro- 
vided by ’’breathing" during thermal cycling, and since 
the cycling frequency and temperature range were dif- 
ferent from previous experience, it was not known whether 
the "breathing" would be adequate to maintain chemical 
stability . 


Present Program 


The present work continues the prior work by means of a 
three-part program: 

Part I - A number of samples of Thermkeep were 
taken from the experimental model for 
chemical analysis to determine whether 
or not the system remained chemically 
stable, and whether or not segregation 


* Registered Trademark Comstock & Wescott, Inc. 
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of components occurred during the prior 
testing program. Thereafter, the Therm- 
keep was melted and partially drained to 
permit visual examination of the upper 
part of the heat exchanger. 

Part II - The Therminol-66 heat transfer fluid 

used in the prior program was replaced 
by Caloria HT^3, and an additional set 
of experimental data obtained for fur- 
ther verification of the computer model 
with a fluid of different physical prop- 
erties . 

Part III - An earlier study had suggested that 

TESS and overall system benefits might 
be obtained in the STETF if a second 
heat exchanger were incorporated in the 
TESS in which the toluene vapor for the 
turbine would be generated, thus elimi- 
nating the external toluene boiler. 
Therefore, a new computer model was 
developed to predict the performance 
of such a TESS. This computer model 
was used to produce a reference design, 
including construction cost, of a unit 
suitable for installation in the Sandia 
STETF. This was compared with the de- 
sign of the single heat exchanger unit 
produced in the prior work. 
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PART I 

ANALYSIS OF THERMKEEP SAMPLES 
AND EXAMINATION OP HEAT EXCHANGER 


The experimental model produced under Contract NAS3'-206l5 
remained in the state in which it was left after several 
weeks of thermal cycling, as described in NASA CR-159^65. 
During this cycling the thermal gradient causes the 
Thermkeep in the upper section of the vessel to be pre- 
dominantly in the liquid state, the Thermkeep at the 
bottom is predominantly solid, with an intermediate con- 
dition between. Since Thermkeep is a non-eutectic 
composition, segregation of components may have occurred. 
In addition, some NaNO^ may have been reduced to NaNO^ 

by the mechanism described below. In order to clarify 
these points, a series of samples of Thermkeep were taken 
for chemical analysis. 


Chemical Analysis of Thermkeep Samples 


The thermal insulation and the external heating shroud 
(used to control the thermal gradient in the insulation) 
were removed from one side of the vessel. With the ves- 
sel at room temperature, seven 10.8 cm diameter holes 
were cut along a vertical line on the vessel wall by 
means of a hole saw. The centers of the holes are 5.1, 
32.7, 60.7, 88.0, ll6.6, 135.5, 170.8 centimeters from 
the bottom of the vessel. These locations span the 
total height of the solid Thermkeep. 

The removal of these sections exposed the solid Thermkeep. 
Approximately 0.23 kg was removed in the form of chips 
through each hole, after which the holes were temporarily 
covered to prevent absorption of atmospheric moisture and 


During the removal of the sample from the third hole from 
the bottom, a very small penetration accidentally occurred 
in one of the heat exchanger coils. However, this was in 
an accessible location and was repaired. This was done by 
flowing nitrogen through the tube while it was being welded 
from the outside. This prevented obstruction of the tube 
bore at the weld. The sampling openings were then sealed 
by steel patches welded onto the tank wall. 
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Examination of Heat Exchanger 


The exposed section of the tank was covered by temporary 
insulation. The heating system which is used for charg- 
ing the TES unit was activated and the Thermkeep was 
completely melted. Approximately 3^0 kg of Thermkeep 
were withdrawn from the tank into a steel drum which 
lowered the level sufficiently so that the upper third 
of the heat exchanger was uncovered. The Thermkeep in 
the drum was covered to prevent absorption of moisture 
and CC> 2 and allowed to self-cool. 

The Thermkeep remaining in the TES unit was allowed to 
solidify and cool. The cover of the vessel was removed 
and the interior of the exposed interior surfaces were 
examined. The heat exchanger was photographed to record 
its appearance. Its condition was compared to photo- 
graphs made of the heat exchanger at the time of manu- 
facture and assembly in the tank. Physically the heat 
exchanger appeared to have undergone no significant dis- 
tortion. The supporting structure appeared to be in 
good condition and was adequate for the planned 
experimental work. The metal surfaces which were ex- 
posed to air above the highest level of the liquid 
Thermkeep were found to have acquired a light coating 
of sodium carbonate crystals. This condition has been 
found consistently in Thermkeep systems on those sur- 
faces which are exposed to ambient air drawn into the 
clearance space by cycling. This material was easily 
brushed off and the underlying steel showed no abnormal 
signs of corrosion. 

As a result of the examination it was concluded that no 
internal repair or maintenance was required before con- 
tinuing with the program. 


Heat Transfer Fluid 


The original program plan was to use Silicone B as the 
heat transfer fluid for the test program. This fluid 
is produced by Dow-Corning and has been designated by 
three different names: Silicone B, Siltherm 800, and 

X2-1162. A quantity of this material was obtained, and 
a qualitative test made to determine whether any obvious 
interaction occurred when Thermkeep and Silicone B are 
in contact at a temperature of 589 K (600 F). It was 
found that when a small quantity of either Thermkeep or 



reagent grade NaOH is dropped into the 589 K (600 F) 
fluid, a vigorous evolution of gas occurs, and the 
temperature of the liquid drops promptly. Analysis 
of the tests indicated that decomposition of the fluid 
was initiated by the NaOH. The decomposition products 
are lower molecular weight cyclic molecules such as 
tetramethyl cyclosiloxane , which are more volatile 
than the Silicone B. These compounds are reported to 
be non-toxic. 

In view of this instability of the Thermkeep-Silicone B 
system, the NASA Project Manager directed that a substi- 
tute heat transfer fluid which is stable in the presence 
of Thermkeep at 589 K (600 F) should be identified and 
substituted, since it is an objective of the program 
not only to verify the computer model of the system 
with a second heat transfer fluid, but also to test 
materials which can be used in a practical system. It 
was subsequently decided that Caloria HT-43 (a product 
of Exxon Co.) should be used in the testing program. 


Chemical Analysis of Thermkeep 


The samples of Thermkeep were analyzed by Skinner & 
Sherman, Inc. for percentage by weight of NaOH, NaNO^, 

NaNC^, MnOg, Na2C0^ , and NaCl. A copy of Skinner & 

Sherman's report is included herewith as Appendix A. 

The values reported are repeated in Table I which also 
shows the percentage of NaOH and the total of the six 
constituents for each sample. 

Samples 1 and 2 were taken from a previously unopened 
drum of Thermkeep flakes and represent the batch of 
material from which the TES unit was originally charged. 
The remaining samples are the pairs which were taken 
from the seven openings in the side of the TES unit. 

All constituents except NaN0 2 were independently deter- 
mined. NaN0 2 was obtained as the difference between a 

determination of total nitrogen (not shown in the report), 
and the NaNO^. The reported results are thus based on 

six independent determinations. 
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TABLE I 


ANALYSIS OF THERMKEEP FROM TES VESSEL, AND ORIGINAL MATERIAL 

(% By weight) 


Sample 

No. 

Location 

No. 

Dist . 
from 
bottom , 
cm 

NaOH 

NaN0 3 

NaNOg 

Mn0 2 

Na 2 C0 3 

NaCl 

Total 

1 

Ong.Mat'l N/A 

8 6.9 b 

8.1+67 

.0306 

.172 

2.62 

2.03 

100.26 

2 



88.74 

7*1+77 

.0307 

.171 

1.1+5 

1.98 

99.85 

3 

1 

172 

91.50 

I+.652 

.0020 

.101 

2.85 

1.73 

100.81+ 

4 



91.55 

5.278 

.0021 

.085 

1.1+8 

1.35 

99-75 

5 

2 

144 

91.17 

5.01+1 

.0023 

.113 

1.55 

1.1+7 

99.35 

6 



91.71+ 

5.326 

.0022 

.110 

1.1+3 

1.63 

100.24 

7 

3 

117 

91.88 

1+.503 

.0019 

.083 

1.38 

1.68 

99.53 

8 



92.23 

1+.377 

.0023 

.110 

1.35 

1.85 

99.92 

9 

k 

89 

89.50 

5.623 

.0023 

.102 

2.89 

I.65 

99.77 

10 



91.66 

I+.878 

.0026 

.095 

1.36 

1.76 

99.76 

11 

5 

62 

94.30 

2.791 

.0017 

.01+2 

1.1+1+ 

1.53 

100.10 

J -C. 

13 

6 

34 

93.93 

2.707 

.0021 

.057 

1.1+0 

1.76 

99.86 

Ik 



9 k . 5 9 

2.5I+6 

.0032 

.053 

1.32 

1.25 

99.76 

15 

7 

6 

89.17 

5.772 

.0098 

.137 

2.88 

1.58 

99.55 

16 



89.88 

5.138 

.0095 

.172 

2.90 

1.1+3 

99.53 



The fact that in every case the deviation of the total 
of the six from 100% was no greater than 0.65% was an 
indication of the precision of the analytical work. 

These totals also indicated that the samples absorbed 
only an insignificant amount of water from the ambient 
air during the sampling procedure. 

As previously described, the samples were removed from 
the TES unit by cutting 10.8 cm diameter circular sections 
from the tank wall. The samples were removed within a 
distance of about 5-1 cm from the tank wall. The samples 
were not mixed, so that the first sample removed may 
represent the material closest to the tank wall, and 
the second sample the material taken from the deeper 
layer . 

Prior to removing the samples, the TES unit was last 
operated in the program described in Report CR-159^65. 

The final tests were simulations of the solar daily 
cycle, the results of which are shown in Figs. 91-100 
of that report. 


Discussion of Analytical Results 


A significant result was that the total sodium nitrate 
content of the samples taken from the TES unit was sig- 
nificantly lower than that in the original material. 

Figures 1, 2, and 3 show the averages of the two samples 
taken from each opening in the tank, plotted against the 
location of the opening. The distance between the cen- 
ters of adjacent openings was 27.6 cm. 

The curves for sodium nitrate, sodium carbonate, manga- 
nese dioxide, and sodium nitrite have similar shapes 
in which the percentage of these constituents tends to 
be lowest at locations 5 and 6 and highest at location 7 
(nearest the bottom of the vessel) . This suggests that 
in a zone about 5.1 cm thick adjacent to the wall of the 
vessel, there was some vertical segregation of these 
components of Thermkeep. 
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Top Location Bottom 

Figure 2. Variation of MnC> 2 with vertical location in tank. 
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Figure 3. Variation of NaNOg with vertical location in tank. 
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In Figure 4 the ratio of the average NaOH content to the 
average NaNO^ content, for each pair of samples, is 

plotted against the sampling location, showing that this 
ratio for all samples is significantly higher than that 
of the original Thermkeep, and that at locations 5 and 
6, it is higher than at the other locations. 

Figure 5 is a phase diagram (Ref. 2) for the system 
NaOH-NaNO^ which shows that the solid which forms upon 

cooling was higher in sodium hydroxide than was the liquid 
phase from which it was deposited. Although this phase 
diagram does not truly represent Thermkeep because of 
presence of 2-3% of NaCl and NagCO^, it is probable that 

during solidification a local variation in composition 
in the vicinity of the heat exchanger coil would result. 
The solid close to the tubes would be expected to be 
lower in sodium nitrate than the solid near the center- 
line of the heat exchanger coils or in the spaces 
between the coils. This might also account for the low 
nitrate content of the Thermkeep which solidified on 
the vessel walls as a result of heat lost through the 
insulation . 


Chemistry of Thermkeep-Steel System 


Studies prior to this contract have revealed some of the 
chemical characteristics of the Thermkeep-steel system 
at 756 K (900 F) . Chemical stability requires contact 
with air. This is provided by the "breathing" of ambient 
air in and out of the clearance space above the Therm- 
keep during thermal cycling due to the changing volume 
of the Thermkeep. The reactions which are believed to 
take place are: 


2 Fe + 6 NaOH + NagFe^ + 3 H ? + 2 Na 2 0 ( 1) 

x Fe + y NaNO^ -* Fe x ° y + y NaN0 2 ( 2) 

v (MnO ) 

y NaN0 2 + | 0 2 - 6 — £> y NaN0 3 (3) 
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2 Na 2 0 + 2 H 2 I •> 4 NaOH 


(4 ) 


The net reaction is: 


( x+2 ) Fe + 2 NaOH + | 0 2 + 2 H 2 0 

(5) 

Na 2 Fe 2°4 + Fe x°y + 3 H 2 


Eq. 1) is the "corrosion reaction." Na 2 Fe 2 0^ is the 

principal corrosion product, and has been identified by- 
X-ray diffraction tests. The presence of hydrogen has 
also been confirmed. Corrosion rates in inches per year 
(ipy) of penetration at 756 K (900 F) have been measured 

at less than 2.5 x 10 ^ cm/year (0.001 ipy). In the pres- 
ent application at about 584 K (591 F), corrosion is 
expected to be insignificant in the time of operation, 
and no measurements were undertaken. 

Eq. 2) represents the formation of iron oxide of uniden- 
tified composition, by reduction of NaNO^. Salt baths 

similar in composition to Thermkeep are widely used for 
metal cleaning, at temperatures of 756-867 K (900-1100 F) . 
These are known to produce oxide coatings on the metal 
by NaNO^ which is reduced to NaNO,-,. 

Eq. 3) represents the oxidation of NaN0 2 to NaNO^ by 
oxygen from the air. This reaction is catalyzed by Mn0 2 . 

Eq. 4) represents the reaction of the Na 2 0 formed by Eq . 

1) with moisture introduced by the "breathing" of humid 
air . 
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Temperature, deg. 


One purpose of the chemical analyses was to determine 
whether sufficient oxygen was absorbed by the system 
to maintain the original level of NaNO^. Since the 

analysis showed that the NaN0 2 in the system was an order 

of magnitude lower than in the original Thermkeep, it was 
concluded that the oxygen absorption rate was sufficient 
to prevent significant conversion of NaNO^ to NaNO,-,. 

The analyses showed that the NaNO^ content of the samples 

was substantially lower than that of the original Therm- 
keep. This raises the question of whether or not nitro- 
gen may have been lost from the system. Based upon prior 
work, it is believed that loss of nitrogen cannot account 
for this phenomenon. 


Prior experiments have been conducted on the Thermkeep 
iron system, at 756 K (900 F), in steel containers which 
were hermetically sealed to prevent entrance of air, in 
which the pressure was monitored. It was observed that 
such a system, at 756 K, initially develops a subatmos- 
pheric pressure which remains stable for about five days, 
after which it rises abruptly. The rise in pressure is 
found to be due to the generation of ammonia. This 
behavior is due to the sequential reduction of NaNO^ — 

first to NaNC^ , then to NH^ , according to the following 
reactions: 


8 Fe + 24 NaOH ->• 4 Na^e^ + 12 H 2 + 8 Na 2 0 


( 6 ) 


3 NaN0 3 + 3 H 2 -*■ 3 NaN0 2 + 3 H 2 0 
3 NaN0 2 + 9 H 2 -> 3 NH 3 + 3 H 2 0 + 3 NaOH 
6 Na 2 0 + 6 H 2 0 -*■ 12 NaOH 


(7) 

( 8 ) 
(9) 
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Summing these equations gives : 


8 Fe + 9 NaOH + 3 NaNO^ -> 

( 10 ) 

4 Na 2 Fe 2 O ij + 3 NH 3 + 2 Na 2 0 


Based upon this chemistry, and the fact that no NH^ has 

ever been detected in the clearance space of a cycling 
system which has had unobstructed access to air, it was 
concluded that nitrogen had not been lost. 

The fact that the NaNO^ content of all samples was sub- 
stantially lower than that of the original Thermkeep 
charged into the TES tank has been attributed to segre- 
gations of chemical components; a small degree of seg- 
regation in the vertical direction and a larger amount 
in the radial direction. 


The radial segregation was believed to be explainable on 
the basis of the NaOH-NaNO~ phase diagram. Figure 5, 
which shows that during the thermal discharge cycles the 
solid which forms on the surfaces of the tubular heat 
exchanger is richer in NaOH than the liquid from which 
it crystallizes, or equivalently the solid has a lower 
percentage of NaNO^ than the liquid. The samples used 

for analysis were removed through openings in the tank 
wall and are therefore solid which formed on or close to 
the wall. Prior to removal of the samples, the last 
operation of the TES unit (Fig. 100, Report NASA CR-159465) 
consisted of a thermal discharge which produced a solid 
layer on the heat exchanger surfaces, leaving liquid in 
the centers of the tubular helices and near the tank wall. 

Upon termination of the test the Thermkeep temperatures 
in the upper two-thirds of the tank were above 565 K, 
which as shown by the phase diagram would contain 50% or 
more of liquid. 
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After completion of the test the circulation of the heat 
transfer liquid was stopped and the TES unit allowed to 
self-cool by losing heat through the insulation. This 
would cause the remaining liquid to solidify starting 
at the tank walls. The solid forming at and near the 
walls would be low in NaNO^, while the remaining liquid 

would increase in NaNO^ until the eutectic composition 

was reached at some distance from the walls. This ex- 
planation is, of course, over simplified because it does 
not take account of the Na^O^ and NaCl which total 

approximately 4% of the original Thermkeep composition. 

A number of tests were considered which might verify the 
presumption that no nitrogen has been lost, and shed 
more light on the segregations. However, the opportu- 
nity to obtain more information on the system in its 
original state was lost when the Thermkeep was completely 
melted for draining and examination of the heat exchanger. 


Summary of Chemical Analyses of Thermkeep 


The following conclusions were drawn as a result of the 
observations described above: 

1. The accuracy of the analyses was good. 

2. There were no significant compounds 
present in the system other than 
those analyzed. 

3. The samples were not significantly 
contaminated with H^O. 

4. No significant NaNO^ was present and, 

therefore, 1) no nitrogen was lost 
from the system, and 2) air "breath- 
ing" is adequate to maintain chemical 
stability . 

5. Variations in the NaOH/NaNO^ ratio 

were a result of segregation of compo- 
nents during cooling. 
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6. A small vertical segregation existed 

as shown by the high NaOH/NaNO„ ratio 
at locations 5 and 6. ^ 

7. Since the NaOH/NaNO^ ratio at all 

seven sampling locations was less than 
that of the original material, there 
must also be a radial segregation from 
the tank wall inward. These radial seg- 
regations had no observed deleterious 
effect on the operation of the TES sys- 
tem in prior tests, and are probably 
acceptable if the system had arrived at 
a stable state as a result of the cycling 
carried out in the prior program of test- 
ing. To determine whether stabilization 
was achieved will require further testing. 

8. If further cycling produces progressive 
segregation to the point where the per- 
formance of the system is unacceptably 
degraded, this condition could be re- 
versed by periodically over-charging the 
TES unit and allowing the liquefield 
Thermkeep to become homogeneous by natu- 
ral circulation of the liquid. 
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PART II 

VERIFICATION OF COMPUTER MODEL 
BY EXPERIMENTAL DATA 


Introduction 


The purpose of this work was to extend the experimental 
verification of the computer model (CM) of the single 
heat exchanger TESS, by means of additional data ob- 
tained from the experimental scale-model. The new data 
were obtained through the substitution of a different 
heat transfer fluid for the Therminol-66 (T-66) used 
previously. The substitute fluid was to have physical 
properties as different from those of T-66 as possible, 
in order to demonstrate that the CM could be used with 
confidence to predict TESS performance for a range of 
fluids . 

The program of testing was to duplicate the prior pro- 
gram, so far as heat transfer rates were concerned, by 
adjusting flow rates to compensate for the different 
physical properties of the new fluid. The fluid tenta- 
tively selected was Silicone B (a product of Dow-Corning) , 
because its physical properties are quite different from 
those of T-66 , and it has superior thermal stability, 
although its cost is substantially higher. When, as de- 
scribed in Part I, it was found to be reactive with mol- 
ten Thermkeep, it was rejected on the basis of safety, 
and Caloria HT-43 was substituted. 


Heat Transfer Fluid 


Caloria HT-43 has a higher specific heat (Cp) but a lower 
density (p) than T-66, and flow rates are inversely pro- 
portional to the product of these values. The following 
compares the two fluids at 316 C (600 F): 


T-66 

HT-43 


P Cp 

g/cm 3 .t/( e .k) P x Cp 

.805 2.63 2.12 

.649 2.93 1-90 
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Since pCp for HT-43 is 90$ that of T -66, the flow rates 
necessary to effect th_e same heat rates at 589 K fluid 
temperatures must be increased by about 10$. The high- 
est steady state flow rate in the test plan was 0.35 
kg/sec (6.6 gpm) for T-66, while the pump is capable 
of 0.40 kg/sec ( 7.6 gpm). The only test in which the 
heat rates were not duplicated was the solar cycle, where 
during charging the HT-43 will reach the maximum pump 
rate sooner than did the T-66. Since the actual flow 
rate is a computer input, this was simulated in the com- 
puter program. 


Experimental Model 


The experimental model is described in detail in NASA 
CR-159465. Briefly, it consists of a steel cylindrical 
thermal storage tank 2.23 m in height and 0.70 m in 
diameter containing 1500 kg of Thermkeep . The heat ex- 
changer consists of 25 helical steel tubes, 0.635 cm 
outside diameter and 23.1 m long, which are manifolded 
at the top and bottom of the tank for parallel flow of 
the heat transfer fluid, downward during charging, and 
upward during discharging. Temperatures are measured 
by 20 thermocouples in wells extending through the 
side of the tank into the Thermkeep to two different 
depths, and additional couples on the outside of the 
tank . 

The temperature of the fluid entering and leaving the 
tank and at other points in the system, as well as 
other system temperatures, are also measured by thermo- 
couples . 

The tank is covered by an inner layer of thermal insula- 
tion 0.457 m thick, and an inner steel shroud carrying 
a grid of electric heating elements. This is covered by 
an additional layer of insulation 0.152 m thick and an 
outer steel shroud. The heated shroud is used to 
establish a desired thermal gradient in the inner insu- 
lation in order to simulate insulation heat loss which 
would occur in a larger unit. 

The auxiliary equipment provides for the supply of heat 
transfer fluid to the TESS at predetermined flow rates 
and temperatures for charging and discharging, at con- 
stant flow rates and at programmed variable rates to 
permit simulation of solar cycles. 
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The experimental model has approximately one-tenth the 
capacity of a TESS suitable for use in the Sandia 
Albuquerque Solar Total Energy Test Facility. The 
height of the model and consequently the length of each 
heat exchanger element are of full scale; the cross- 
section and consequently the number of heat exchanger 
elements are one-tenth of full scale. 


Preparation of Test Unit 


After samples were taken from the storage unit, about 
430 kg (950 lbs) of Thermkeep flakes were added to re- 
place what had been drained to allow visual inspection 
of the heat exchanger. During previous testing the 
liquid level of Thermkeep was about 15 cm lower than 
necessary to cover the entire heat exchanger at the 
maximum temperature, since a large portion of Thermkeep 
usually remains solid at the bottom during normal duty 
cycles. Therefore, enough Thermkeep was added to 
increase the average level by 15 cm. 

The storage tank was closed and covered with insulation. 
The heated shroud was then replaced and its heating 
elements checked and reconnected. The purpose of the 
heated shroud was to produce an elevated "ambient" 
temperature for the tank. The shroud was turned on and 
over a period of ten days was balanced so that the top, 
sides, and bottom were at an average temperature of 
about 493 K. 

All of the original heat transfer fluid (T-66) was 
drained from the system and air was blown through the 
pipes and coils for 24 hours. The system was then 
charged with 10 gallons of Caloria HT-43 (HT-43) which 
was pumped through all of the plumbing, collected, and 
circulated once again. Then the HT-43 was drained and 
air was pumped through the pipes for 24 hours. The pro- 
cedure was repeated with fresh HT-43. After draining 
the T-66 and after each HT-43 flushing, the water- 
immersed cooling coil was removed, inverted, and 
drained. This coil was the only location where more 
than one pint of liquid could collect. The liquid con- 
tent of the heat exchanger was measured between flush- 
ings by balancing its pressure head against a water manom- 
eter, on the upper and lower manifolds, and found to be 
no more than 250 cmh Therefore, in the worst case, 
each 38 liter (10 gal) flush was at least a 100 to 1 
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dilution, and with the final charge of 265 liters (70 
gal)., another 1000 to 1 dilution. 

The unit was heated to various states of charge and the 
liquid level of Thermkeep was measured with a dipstick. 
This showed that the unit contained 1515 kg of Thermkeep, 
or about 75 kg (165 lhs) more than before disassembly, 
and that none of the heat exchanger surface was exposed 
above the Thermkeep level during normal cycling. 

Two minor problems were encountered before the start of 
tests: The cooling water barrel had several leaks 

and was replaced, and thermal siphoning to the expansion 
tank became an increasing problem. The siphon pathway 
was traced to a broken relief valve between the heating 
reservoir and the expansion tank, and the valve was 
removed. An alternate escape pathway already exists 
with a hand-operated valve, and this was fixed in the 
open position. 

All thermocouples were found to be in working order and 
were connected to a 60-channel "Datalogger" (FLUKE Model 
2240B), which automatically records temperatures on 
paper tape, on a predetermined schedule. 


Testing Procedure 


The experimental data taken include the flow rate of the 
heat transfer fluid, HT-i|3, its temperature entering and 
leaving the TES unit, and temperatures of the Thermkeep 
in the TES unit measured by sets of thermocouples in 
wells extending into the tank, in two vertical rows, one 
to a depth of 11.0 cm and the other to a depth of 27-3 
cm from the tank wall (the tank diameter is 69.9 cm in- 
side) . 

The temperature data, recorded on paper tape by the Data- 
logger, were manually transferred to files of a small 
computer. The computer model (CM) of the TESS required 
a larger computer located in an outside facility for the 
CM performance predictions. The results of these predic- 
tions were also transferred to the small computer file, 
which was then used to print graphical displays of the 
experimental data and superposed CM predictions. In 
these graphs, the maximum deviation of the printed 
point from the true point is 2.5°C. 
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Inputs into the CM program include the physical proper- 
ties of the heat transfer fluid, and its flow rate. 
Twenty- five elements are used in the analysis. The 
initial condition data for the CM includes the tempera- 
ture of the Thermkeep in each element, as measured at 
the start of the corresponding test. 

The experimental data and the CM predictions of perform- 
ance of Thermkeep temperatures vs. time are shown in 
several sets of computer-printed graphs which are to 
be found in Appendices B and C. Temperatures were re- 
corded and graphed at half-hour intervals, but the data 
in Appendices B and C are abbreviated to conserve space. 
Superposed on the data are the temperatures predicted 
by the computer model. The graphs show the Thermkeep 
temperatures vs. the location in the tank, measured 
from the top of the tank. Also shown is the location 
of the upper surface of the Thermkeep as calculated by 
the CM; this location depends upon the temperature vari- 
ation within the tank and the fraction of Thermkeep 
which is solid. 


Enthalpy-Temperature of Thermkeep 


The CM uses a table relating the enthalpy and tempera- 
ture of the Thermkeep, which is used to determine the 
energy of each of the finite elements. Data for the 
table were to have been determined under a subcontract 
of the preceding contract. This failed to produce use- 
ful information and in its place data obtained by 
Comstock & Wescott during large scale cooling tests 
were used. These data were input to the CM in the form 
of a table. The table used in the prior work to corre- 
late data and CM predictions is identified as "Table 4" 
and is shown graphically in Figure 6. 

The ability to produce large numbers of graphical dis- 
plays of superposed experimental data and CM predictions 
quickly and inexpensively by computer data processing 
enabled a more complete analysis of the temperature- 
enthalpy relation than was possible during the preced- 
ing study. This showed that Table 4 is a less than ideal 
temperature-enthalpy representation, and also that some 
degree of hysteresis occurs in the latent heat range 
during heating-cooling cycles. In order to study these 
effects, the enthalpy table was adjusted by trial-and- 
error until two new tables were produced, one which 
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Enthalpy, kJ/kg 





correlates the CM predictions better with cooling data, 
and another with heating data. These are shown graphi- 
cally in Figure 6, where Table 7 is used for heating, 
and Table 9 for cooling. 

In addition. Table 10 was constructed and used as a 
compromise to correlate both heating and cooling data. 


Hysteresis Study 


The construction of these new tables was based on CM 
predictions correlated with experimental data obtained 
with a constant flow rate of the heat transfer fluid 
HT-43 at 0.125 kg/sec (3 gpm). Thereafter, correla- 
tions at 0.042 kg/sec (1 gpm) for heating and cooling 
and at 0.305 kg/sec (7-3 gpm) for cooling were produced. 

Table II shows the conditions of nine of the tests which 
were run, each of which ran for a total of six hours, 
except for No. 9 which ran four hours. Data were taken 
at the start of each test and thereafter at one-half 
hour intervals. Correlations for each data set were 
made with CM predictions. From these, the set showing 
the greatest deviation, and the starting condition, were 
selected. The graphs for these are to be found in 
Appendix B. In these graphs the symbols have the fol- 
lowing meanings : 

A = temperatures measured In deep wells, 

o = temperatures measured in shallow wells, 

+ = temperatures prediced by the computer 
model . 


Figures B-l - B-4 show the deviation using Enthalpy Table 
4. The effect of hysteresis (which is not accounted for 
in the CM) is shown by the fact that Thermkeep tempera- 
tures in the temperature range of latent heat are higher 
than predicted during heating, and lower during cooling. 

Figures B-5 and B-6 show considerably better correlation 
for heating through the use of Enthalpy Table 7 in the 
CM, rather than Table 4. 
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TABLE II 


TEST CONDITIONS FOR HYSTERESIS STUDY CORRELATION 


Test 

No. 

Figures of 
Appendix B 

Mode 

HT-4 3 

Flow Rate 

Enthalpy- 

Table 

1 

B-l, B-2 

Charge 

0.125 kg/sec 
(3 gal/min) 

4 

2 

B-3, B-4 

Discharge 

0.125 kg/sec 
(3 gal/min) 

4 

3 

B-5, B-6 

Charge 

0.125 kg/sec 
(3 gal/min) 

7 

4 

B-7, B-8 

Discharge 

0.125 kg/sec 
(3 gal/min) 

9 

5 

B-9, B-10 

Charge 

0.042 kg/sec 
(1 gal/min) 

10 

6 

B-ll, B-12 

Discharge 

0.042 kg/sec 
(1 gal/min) 

10 

7 

B-13, B-14 

Charge 

0.125 kg/sec 
(3 gal/min) 

10 

8 

B-15, B-16 

Discharge 

0.125 kg/sec 
(3 gal/min) 

10 

9 

B-17, B-18 

Discharge 

0.305 kg/sec 
(7.3 gal/min) 

10 
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Figures B-7 and B-8 show better correlation for cooling 
through the use of Enthalpy Table 9 in the CM, rather 
than Table 4. 

Figures B-13 - B-l6 show the maximum deviation obtained 
with Enthalpy Table 10 used both for heating and cool- 
ing. 

The remaining figures show the maximum deviations meas- 
ured at higher and lower flow rates, using Enthalpy 
Table No. 10. 

On the basis of these tests, it appears that the use of 
Table 10 for both heating and cooling provides an ade- 
quate representation of performance for the purposes of 
this program. Table 10 was subsequently used to cor- 
relate all data from the heating and cooling cycling 
tests, and selected data obtained during the prior pro- 
gram with T-66 as the heat transfer fluid. 

Some preliminary consideration has been given to the 
problems of incorporating into the computer model a 
means for accounting for the apparent hysteresis effect. 
However, this appears to be very complex, and to require 
quantitative information not now available. The true 
nature of the phenomenon is not understood. A super- 
cooling and/or superheating may be involved in which the 
temperatures of latent heat effects are elevated during 
heating, and/or depressed during cooling. The degree 
of elevation or depression may be functions of the 
temperature and/or the rate of heating or cooling. A 
contributing factor may stem from the fact that upon 
cooling from the liquidus range, the first solid which 
forms on the heat exchanger is pure NaOH whose melting 
point is given by literature references as 593 K 
(cf Fig. 5), which is above the temperature of the heat 
transfer fluid entering during charging. Upon reheat- 
ing with a heat transfer fluid whose temperature is 
several degrees below its melting point, reliquefying 
must be principally by solution in the liquid phase at 
the liquid-solid interface rather than at the heat ex- 
changer surface. (The CM is based upon this mechanism.) 

Experimentally determined information would be needed 
relating the rate of change of the enthalpy of the 
Thermkeep to its temperature, its rate of heating or 
cooling, and, if necessary, its immediately preceding 
temperature history. This information may be extract- 
able from the experimental data already in hand, which 
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includes an almost continuous temperature history of 
Thermkeep temperatures throughout the TES unit for a 
wide range of heating and cooling rates, and for cycles 
in which heating is immediately followed by cooling, 
and vice versa. 

If the necessary information can be obtained, it must 
then be introduced into the CM in a way which follows 
the first law, the conservation of energy. That is, 
when an element of Thermkeep is carried through a cycle 
within the temperature range of the hysteresis effect, 
the cycle must form a closed loop on the temperature- 
enthalpy diagram. This cannot be accomplished simply 
by using two tables, one for heating and one for cool- 
ing and shifting from one to the other. For example, 
with reference to Figure 6, if Table 7 were used for 
heating and Table 9 for cooling, and if an element were 
heated from 558 K to 573 K along curve 7> then cooled 
along curve 9 back to 558 K, more energy would be de- 
livered by the element during cooling than would be 
received during heating. Upon repeated cycling between 
these temperatures, the element would progressively 
lose energy. 

This might be avoided by introducing a third t-h rela- 
tionship such as Table 10 which represents a zero-rate 
relation between t and h, and by requiring that the 
maximum and minimum temperatures of cycles which are 
within the hysteresis temperature range must lie on this 
curve . 

However, correlations next to be described, in which Enthalpy 
Table 10 was used for all cycles are believed to be ade- 
quate for the intended uses of the CM. 


Data Correlation 


The data correlations consisted of comparisons between 
experimental data (DATA) and computer predictions (CMP). 
The correlations include Thermkeep (TK) temperatures 
and the temperatures of the heat transfer fluids enter- 
ing or leaving the upper end of the heat exchanger (TU) 
and entering or leaving the lower end (TL). 
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Correlations were made using data for Therminol-66 
(T -66) from the prior program and for Caloria HT-43 
(HT-43) obtained in the present program. The results 
are given in Appendix C in graphical form. 

Figures C-l through C-12 show HT-43 inlet and outlet 
temperatures vs. time. Figures C-13 through C-48 
show the corresponding Thermkeep temperature profiles 
in the tank. 

Figures C-49 through C-59 show T-66 inlet and outlet 
temperatures vs. time. Figures C-60 through C-75 show 
the corresponding Thermkeep temperature profiles in 
the tank. 

The tests performed consisted of four types: 

1. Single discharge, in which the TES is 
initially charged, and is discharged 
by a constant flow of fluid. 

2. Single charge, in which the TES is ini- 
tially in a partially charged state, as 
at the end of a typical cycle, and is 
charged by a constant flow of fluid. 

3. Cyclic, in which a cycle consists of a 
4-hour charge and a 4-hour discharge at 
constant fluid rate. 

4. Solar cycle, which conforms to the solar 
daily cycle described in NASA CR-159465. 

The solar daily cycle is shown in Fig- 
ure 7 • 

The conditions of the tests performed in the present pro- 
gram are given on Table III. 

The conditions of the prior tests with T-66 are shown on 
Table IV. The correlation of the data with CMP is given 
in NASA CR-159465, which used the original Enthalpy Table 4. 

Correlation of the original data for T-66 temperatures 
with CMP was repeated for inclusion herein, using Enthalpy 
Table 10, and also TK temperature correlations for Test 
018 and Test 023. 
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Figure 7. Typical solar daily cycle at Sandia Total Energy Test Facility. 



TABLE III 


TESTS WITH CALORIA HT-43 
UNDER PRESENT PROGRAM 


Test 

No. 

Type of Test 

HT-43 Flow 
kg/sec 

Rate 

gpm 

201 

Single discharge, 6 hrs. 

0.14 

3-3 

202 

Single charge, 6 hrs. 

0.14 

3.3 

204 

Single charge, 6 hrs. 

0 .16 

3.7 

214-2 

Single discharge, 6 hrs. 

0.046 

1.1 

215-2 

Single charge, 6 hrs. 

0.046 

1.1 

218 

Cyclic; 4-8 hr. cycles 

0.080 

1.9 

219 

Single discharge, 4 hrs. 

0.32 

7-7 

223 

Solar cycle; 2-12 hr. cycles 

Variable 


TABLE IV 

TESTS WITH THERMINOL-66 
UNDER PRIOR CONTRACT* 


Test 


T-66 Flow 

Rate 

No. 

Type of Test 

kg/sec 

gpm 

001 

Single discharge, 6 hrs. 

0.15 

2.9 

002 

Single charge, 9 hrs. 

0.15 

r\j 

CO 

004 

Single charge, 6 hrs. 

CO 

i — 1 

o 

3.4 

014 

Single discharge, 6 hrs. 

0.053 

1.0 

015 

Single charge, 6 hrs. 

0.053 

1.0 

018 

Cyclic; 3-8 hr. cycles 

0.106 

2.0 

019 

Single discharge, 3 hrs. 

0.35 

6 . 6 

023 

Solar cycle; 2-12 hr. cycles 

Variable 

* See 

NASA CR- 159 465 
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For the constant rate charging and discharging tests 
with HT-43 3 TK temperature profiles are shown for the 
initial state, and after 3 hours and 6 hours, except 
in the case of Test 219 which was of shorter duration. 
After 6 hours, correlation is generally good for flow 
rate of 0.14 kg/sec. An exception is noted in the 
case of Test 204, at 0.16 kg/sec, which appeared be- 
tween the 4th and 5th hours, the cause of which is 
unexplained. The correlation of HT-43 temperatures 
for all these tests is good, with a maximum discrepancy 
t the end of the discharging tests of 6 K (Test 201). 
The HT-43 temperature correlation of Test 204 is good 
(4 K maximum discrepancy) in spite of the discrepancy 
noted above in the TK temperature correlation. 

The high constant discharge rate test (Test 219) showed 
more discrepancy in the TK temperatures than the lower 
rate tests. However, The HT-43 temperature correlation 
is exceptionally good, the maximum discrepancy being 
only about 1 K. 

The cyclic test with HT-43 (Test 218) consists of 4 
cycles, each consisting of 4 hours of discharging fol- 
lowed by 4 hours of charging at a constant flow rate 
of 0.080 kg/sec (1.9 gpm) . During the discharge por- 
tions of this test, the same type of discrepancy occur- 
red as in Test 204. This is reflected in the HT-43 
temperature graphs; however, the maximum discrepancies 
between DATA and CMP were less than 10 K at the end of 
32 hours of cycling. 

Test 223 shows the results of the solar cycle test with 
HT-43. Although some discrepancy is apparent in TK 
temperatures between 6 and 8 hours, they disappear at 
10 hours. Deviations between DATA and CMP for HT-43 
temperatures are no greater than about 5 K throughout 
the 36 hour duration of the test. 

The TK temperatures of the cyclic Test 018 using T-66 
show better correlation than that of Test 218 using 
HT-43, and this is reflected in a very good correla- 
tion for the Test 018 fluid temperatures. On the other 
hand, the correlation of fluid temperatures for the 
Solar Cycle tests (Tests 223 and 023) is better for 
HT-43 than for T-66. The Solar Cycle test with T-66 
shows approximately the same degree of correlation 
with Enthalpy Table 10, as with the previously reported 
Enthalpy Table 4. 
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In view of the generally good correlation of the fluid 
temperatures of the steady flow rate tests, it is sur- 
prising that the discrepancy is so high in the solar 
cycle with T-66. It is noted that the large discrepancy 
occurs during the decreasing fluid flow portion of the 
charge phase of both the first and second cycles, and 
that the correlation is good during the discharge por- 
tion of both cycles (hours 10-12 and 34-36) which fol- 
lows the charge portions (hours 3-9 and 27-33) • Since 
the fluid flow rate during charging is variable in the 
solar cycle (starting and ending at substantially zero), 
it may be that the experimental model does not perform 
as assumed in the computer model. The latter assumes 
uniform flow through all 25 heat exchanger elements and 
that there are no lateral temperature variations within 
the elements. There is evidence that this was not the 
case during the periods in which the discrepancy occurs. 

Figure 75 of Appendix C shows temperatures at the surface 
of the tank (□ ), in the shallow wells (o), and in the 
deep wells (A), at 32 hours of the solar cycle Test 023, 
using T-66. Temperature variations of as much as 25 K 
are noted between the surface and deep well locations, 
and about 13 K between the deep and shallow well loca- 
tions. Similar variations were recorded at the 9th hour 
of the test which is the corresponding point in the first 
cycle. A similar pattern of variations is apparent in 
the solar cycle test with HT-43. 

These temperature variations may result from variations 
in fluid flow rates through the corresponding elements 
of the heat exchanger. If this is the case, it will un- 
doubtedly be correctable by manifolding design changes. 


Summary 


The practical significance of this work lies in the 
correlation of the measured heat transfer fluid outlet 
temperature with that predicted by the computer model. 
For the six constant rate tests of charging and dis- 
charging with HT-43, the average difference between 
experimental and predicted temperatures is 3 K, and the 
maximum is 6 K. For T-66, the average is 2 K and the 
maximum is 7 K. For the constant flow cycle tests 
(218 and 018 ), the discrepancy at the end of 3 cycles 
is 6 K for HT-43 and 2 K for T-66. For the solar cycle 
with T-66, a maximum discrepancy of 14 K was noted, 
while for HT-43, the maximum was only 5 K. 
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It is believed that the computer model is reliable enough 
to be used in the design of phase-change TES systems, 
using Thermkeep as the heat storage medium, and non- 
phase-change heat transfer liquids which are generally 
similar in properties to HT-43 and T-66 . 
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PART III 

DUAL HEAT EXCHANGE DESIGN ANALYSIS 


Introduction 


Early investigations into the use of Thermkeep* as a 
phase-change heat storage medium for solar power genera- 
tion systems suggested that a more economical storage 
unit could he produced if both the collector fluid (T -66) 
and the working fluid (Toluene) of the power cycle were 
caused to flow through the storage unit. This would 
eliminate the cost of a separate boiler but would in- 
crease the cost of the heat exchanger within the stor- 
age unit . 

It was felt that the thermodynamics of the storage and 
retrieval process would be superior in the case of the 
two-fluid storage unit because of the inherent irreversi- 
bilities and large AT’s associated with the fluid pinch- 
points which occur when heat is exchanged between the 
non-phase change fluid (T-66) and the phase-change stor- 
age medium and power fluid. Early attempts to quantify 
this did indeed indicate a substantial benefit in terms 
of overall storage unit size for acceptable perform- 
ance with the use of a two-fluid storage device. 

For this program, an analysis was developed which account 
for the use of two fluids in a typical daily duty cycle 
according to which the solar collector fluid flows dur- 
ing hours of collection at a rate determined by the inso- 
lation intensity. Moreover, the power fluid, which is 
a vaporizing fluid such as Toluene, flows during the 
hours of power generation. Consequently, at some times 
both fluids flow; at other times only one flows. 


Approach 


There are two ways to approach the analysis and the first 
effort in the formulation assessed the relative complica- 
tions and benefits of the two approaches. As might be 
expected, the approach which is more precise is consid- 
erably more difficult to perform. 


* Registered Trademark Comstock & Wescott, Inc. 



The simpler approach is to view the heat exchange pro- 
cess as one in which each fluid exchanges heat with 
the storage medium, but not directly with the other. 

A heat balance is made on the numerical element of 
storage medium at the end of a time step during which 
the temperature of the medium is assumed to remain 
constant. Then an updated temperature distribution 
is produced. 

The other approach which requires additional simpli- 
fying assumptions in order to be made solvable is 
to consider also the direct heat exchange between 
the two fluids. This permits a substantial amount 
of heat transfer from the collector fluid to the 
power fluid without the intervening resistance of the 
storage medium. It also reduces the quantity of 
medium required because the medium need have capacity 
only for the heat stored, rather than all the heat 
flowing from the collectors to the power fluid. 

Direct heat exchange between the collector fluid and 
the power fluid, when both fluids are in energy bal- 
ance, has a clear thermodynamic advantage, albeit one 
which has not yet been quantified. In the existing 
system configuration, whereby the boiler and storage 
system exist as two separate components, collector 
fluid heats the power fluid directly during a substan- 
tial part of the daily cycle. Since all the fluid flow 
in the system under analysis herein must flow through 
the storage unit, it seemed that a direct heat exchange 
path between the two fluids would be preferred over 
one which imposed the heat transfer resistance of the 
storage medium. 

Therefore, the analysis has been formulated to account 
for simultaneous heat transfer between the fluids and 
between each fluid and the storage medium, i.e., three- 
way heat exchange. The purpose of a parametric study 
would be to evaluate the effect of changes in the heat 
transfer path resistances between the three media. More- 
over, the results of optimization in this configuration 
would be compared in terms of cost, performance, and 
operational characteristics to the current standard de- 
sign which passes only collector fluid through the 
storage unit and uses a separate boiler to heat the 
power fluid. 
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General Features 


The prior analysis may be referred to as a design analy- 
sis — that is to say it Includes a rather detailed 
series of calculations relating to a number of specific 
characteristics of the storage system. It works from 
an input set which consists of basic geometry from which 
it derives heat transfer parameters such as surface 
areas and pressure drops. It does a cost analysis. It 
includes axial conduction heat flows and losses to the 
environment. It is limited, however, to a specific 
geometry — that of coiled tubes in a parallel configu- 
ration . 

At the outset of the present work, it was not known 
what sort of geometry would be preferred for a three- 
way system. A determination was made to provide a less 
detailed analysis which would describe the basic heat 
transfer and thermodynamic behaviour of the process and 
to overlook for the moment other effects which have al- 
ready been shown to be relatively insignificant. On 
this basis, calculations with regard to geometries, 
pressure drops, losses to the environment, and axial 
losses have been put aside. 

It will not be difficult to introduce these later if 
wanted. It will require, however, a clearer, more 
limiting description of the desired design configura- 
tion, one not deemed appropriate at this time. Also, 
it was not entirely clear at the outset that a problem 
of such complexity as this would lend itself to a use- 
ful solution, and introduction of additional details 
was not considered prudent. 


General Description 


The analysis considers variations in the temperatures of 
the three media to be one-dimensional; temperatures vary 
only in the direction of flow of the fluid media. The 
hot end of the unit is the top, as shown on Figure 8. 
This stabilizes the thermal gradient established in the 
storage medium by preventing mixing due to natural con- 
vection flow. 

Heated fluid from the solar collector, assumed to be a 
liquid with a constant specific heat, C^, flows through 
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Transfer fluid pump Power fluid pump 


Figure 8. Solar electric power generating system schematic. 
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the unit from top to bottom. Its inlet temperature, 

assumed to be constant, is T t . while its exit 

y in 

temperature is T L ^ . Its flow rate, also a time 
variable, is W p , its value being a function of the col- 
lector heat output and the fluid temperature. 


Cooled power fluid, assumed to be a vaporizing liquid 
with constant specific heats in the sensible regions, 
flows into the unit at the bottom and out at the top. 
Its liquid phase specific heat is C^, its vapor phase 

specific heat is C^., its heat of vaporisation is Q^, 

its flow rate is Wp (assumed to be constant in time), 

its inlet temperature is Tp (constant in time) and 

its exit temperature is Tp t which is a time-varying 

output of the analysis. 

The power conversion system, however, in order to pro- 
duce a constant output, P, requires a constant inlet 
temperature, T p des , and a boost heater is assumed to 

make up the difference between the desired and the 
actual delivery temperatures. In general, this supple- 
mental heat is 


Q supp W P C V ^ T P,des 


T P,out ) * 


The solar collectors produce a scheduled amount of heat, 
Q co ii> which is assumed to be a function of time of day 

but repeatable from day to day. The collector fluid 
attempts to absorb this heat by varying its flow rate 
according to 



Q 


coll 


C L (T L,in 


L,out ) 
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At some point, the fluid pump reaches a maximum flow, 

W T , and any heat which cannot be absorbed due to 
L,max 5 J 

this limiation is assumed to be lost through a heat re- 
jection device. It is calculated by 


Q lost Q coll W L,max C L (T L,in T L,out^ 


Both Q lost and Q SU p p are important operational charac- 
teristics of the system. 


The task 
eters Q 


of 


the 

P, 


erties 


T. 


coll 5 
and the 
T, 


L,out’ P,out 5 


analysis is 

T L,in> T 

storage 

Q-, , Q 

lost, supp. 


to use as input the 
the fluid thermal 


P , in 5 

unit geometry, and to 
and the storage 


param- 
prop- 

calculate 

medium 


temperature variation in space and time. The material 
discussed so far, the parameters external to the stor- 
age unit itself, constitutes the essence of the MAIN 
portion of the computer analysis. The major portion of 
the computer analysis, a subroutine called ELANAL, deals 
with the processes within the storage unit. It receives 
the flows and inlet temperatures of the fluids from the 
MAIN program. 


Heat Exchange Analysis 


Description 


The analysis of the heat exchange processes within 
the storage unit rests upon the same assumptions 
made to develop the heat exchange analysis for the 
single fluid case. The actual procedure is greatly 
complicated by the fact that two fluids flow counter 
current and that one fluid exists in three states: 
subcooled liquid, two-phase boiling fluid, and super 
heated vapor. This last factor imposes the defini- 
tion of three distinct spatial regions: Region I 

where the power fluid is superheated. Region II 
where boiling is taking place, and Region III where 
the power fluid is subcooled. 



Two basic approximations are made in the solu- 
tion methodology; finite differencing in time, 
and finite elements in space. With regard to 
the temporal variation in the parameters, the 
approximation is made that the time axis is 
divided into sufficiently small time steps that 
the system may be viewed as stationary during 
said time steps. 

Heat flows are computed from temperature poten- 
tials assumed not to change significantly during 
the time step. The heat flows contribute to a 
stepwise change in the thermal state of the 
storage medium, the result of which is assumed 
staionary for the next time step. The storage 
medium is assumed to have the only thermal 
capacitance of any significance. 

Thus, the differential equation 


Q = MC 


dT 

dt 


is approximated by 


AT 


_Q 

MC 


At 


where the AT of a piece of storage medium of mass, 
M, and specific heat, C, is that created by impos- 
ing a fixed heat flow Q on it over a time step. 

At. This is the finite differencing in time. 

Spacial variation of the temperature of the stor- 
age medium is approximated by a stepwise variation 
in the (vertical) X-direction. The unit is divided 
in the X-direction into N £ elements of size AX and 

of mass AM. Figure 9 describes an element. Each 
element is assumed to be of uniform temperature, 

T^. For each element, heat flows among the media 

are computed from the knowledge of the inlet tem- 
peratures of each fluid, T n . . , and T 0 the 

ri , 1 + 1 0,1 
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element temperature, and the heat flow re- 

sistances between the media. 

The heat flow between the fluid media, ^ , 

contributes a portion of the heat needed to 
change the fluid temperatures (or states). The 
heat flows between the storage medium and each 
fluid contribute the remainder. These heat 
flows, Q he i and Q EC ^ , also, when summed to- 
gether constitute the net heat flow to the 
storage medium finite element and, when multi- 
plied by the appropriate time step, t, determine 
the element temperature change. Consequently, 
the updated set of element temperatures may be 
calculated to be used in the next time step. 

The three-way heat flow analysis for the power 
fluid is different in the sensible regions 3 where 
temperature changes with heat exchange, from that 
in the boiling region, where quality changes with 
heat exchange. Consequently, the normally equally 
sized elements must be divided into fractional 
elements when a regional transition occurs, 
each fraction being described by its appropi- 
ate relationships. 


Governing Equations 


For any element with constant element tempera- 
ture and fluids flowing counter-current ly , the 
following equations describe the heat transfer 
and thermodynamics : 


dT H (X) 


dQ RE + dQ 
W L C L 


HC 


dT c (X) 


dQ 


EC 

w p c c 


dQ 


HC 


H7 



dQ 


dQ 


dQ 


HE 


HC 


EC 


(UA) he (T r - T)dX 

Cua) hc ct h - T c )dX 

(UA) ec (T - T c )dX 


where (_UA ) EE is the product of U, the overall 

heat transfer coefficient, and A, the effective 
surface area for heat transfer, associated with 
the heat flow Q^ E and likewise for (UA)^ and 

( UA ) E q . The UA are derived in the normal way con- 
sidering the overall heat flow path between the 
two relevant media. 

There are two solution forms of use to this situ- 
ation. In the general case the power fluid is in 
the sensible region and is, in fact, a function 

of X. In this case, the solution is of the form 


T C (X) 



(1) 


t h (x) 




( 2 ) 


where A^, A^ are solutions to a characteristic 
equation; L is the effective element length; , 

Zg, Z^j and Z^ are groups of terms formed from the 
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UA, W L C L , WpC^ , and the element temperature, T, 

which is constant. and C ^ are constants 

determined by substitution of the boundary con- 
ditions, i.e., the known inlet temperatures. 

For the more specific case of T (X) being con- 

stant, which describes the boiling nrocess 
occurring at saturation temperature, the solu- 
tion for T^(X) is of form: 


■'X^X 

T h (X) = Z^e\ L /+ Zg. 


Heat flows for all paths may be calculated by 
integration of the heat flow equations once the 
temperature expressions are known. For example. 


Q HC ^ ua ^hc 



(T h (X) - T c (X))dX 


Application 


In subroutine ELANAL, the element heat flows are 

determined from relationships arising from the 

basic equations discussed above. It is to be 

noted, however, that the fluid temperatures enter 

ing each element are initially unknown. In fact, 

all that are known are temperatures, T p , = T . 

o j i p 3 in 

and T tt = T t . from the external system. 

Not only that, there are three distinct regions, 
previously described, the boundaries of which are 
not known until solution is achieved. 



It can be shown that a set of elements totally 
composed of fluids in uniquely sensible regions 
can be solved for directly through simultaneous 
solution of a set of linear equations in the 
unknown fluid temperatures at the element 
boundaries . 

This is suggested by the fact that elements 

yield 2N E equations for temperatures of which 

there are 2N fi + 2, of which 2 are the 

known inlet temperatures for each fluid. 

In fact, the matrix of coefficients, if the 
equations are properly arranged, is a tridiagonal 
matrix, and special computer subroutine packages 
exist for this situation. 

The computer analysis proceeds in the following 
way : 

. . A position is assumed for the point where 
the power fluid reaches a quality of 1.0. 

Here, the temperature is T sai -j the satura- 
tion temperature, and all elements above 
this point lie in a region where the sen- 
sible heat equations govern. The inlet 
temperature to this set of elements is 
known for the collector fluid as well since 


. . The aforementioned solution is used for the 
set of linear equations. This yields all 
the desired fluid temperatures including 
the collector fluid temperature at the X- 
coordinate corresponding to power fluid 
quality of 1.0, hereafter called T^p. 


. . Knowing the two associated fluid temperatures 
at the end of the boiling region, the heat 
flow equations for T c = constant may be used 

to work backwards to the point where quality 
is exactly 0.0. Neither this point nor the 
point where quality is 1.0 necessarily fall 
on an element boundary. Not only will the 
point be determined but also the associated 
collector fluid temperature, called hereafter 
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. . The region below this point will consist 

of power fluid entirely in the sensible 

region and may be solved for again by 

the set of simultaneous linear equations 

and knowledge of T t T and T D . . This will, 

J_J-L L y in 

however, result in a calculated temperature 
for the power fluid at the outlet, where a 
quality of 0.0 was expected, which will dif- 
fer from T sa £. The extent to which it dif- 
fers is a measure of the accuracy of the 
choice of the point where quality is iden- 
tically 1.0. 

There are two approaches to the use of this method. 
One may be called a closed loop approach and the 
other an open loop approach. The closed loop 
approach is an iterative approach whereby knowledge 
of the error in temperature match where quality is 
to be 0.0 is used to adjust the point where quality 
is 1.0. Adjustment proceeds until the error falls 
within acceptable limits. This iterative approach 
is more accurate but more costly in terms of compu- 
ter processing time. 

The open loop approach uses a method of prediction 
of the point where quality is 1.0 and assumes that 
by proper selection of the prediction method a 
satisfactory accuracy will be achieved. The cost 
of such an approach is much less but accuracy is 
apt to be sacrificed. This is the approach taken 
and results to date are considered satisfactory. 

A method of tracing the movement of the saturated 
vapor point has been developed which combines a 
knowledge of the current error with a predictor of 
the expected result of having updated element 
temperatures . 

Appendix E is a listing of the Fortran program used 
for this work. 


Performance Analysis 


A design is specified by certain limited input parameters. 
The surface areas associated with each heat flow path, as 
well as heat transfer film coefficients for each fluid 
define the geometry. The collector fluid film coefficient 
is input as a function subprogram whose value depends upon 
fluid flow rate. Since the power fluid flow rate is fixed 
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over a cycle, the only variation provided by the analy- 
sis is a dependence upon quality in the two-phase region. 

The duty cycle is, as in the single fluid situation, 
such that the solar collectors produce a normal daily 
variation in output, peaking at solar noon. The power 
system operates for a certain part of the cycle, in the 
early morning and early evening, when the collector out- 
put is zero. Overall the collectors can provide more 
over a day than the power cycle needs and storage smooths 
the daily variation in collector output in an attempt to 
provide constant heat input required. 

Each design run is accomplished by setting the unit at 
full charge, i.e., uniformly at peak cycle temperature, 
and letting the daily cycles progress until heat input 
and output are equal within acceptable limits. This 
stable repetitive operation is the mode of operation 
assumed to be characteristic of the design. Once sta- 
bility is achieved, then the important rating criteria 
can be established. 

When equilibrium is established, the capacity of the unit 
is available as the daily "Storage Output." A "Figure 
of Merit" is defined as the Storage Output divided by the 
sum of the Storage Output and the supplemental heat , 

which is defined above. The "Specific Cost" of 

the TES unit is defined as the total cost divided by that 
portion of the daily output which passes through the stor- 
age medium. 

An analysis has been used to determine the expected per- 
formance of such a system and the results appear in the 
sections to follow. Also, the optimal design from earlier 
work for a single fluid system has been rerun with the 
higher order effects such as insulation heat loss and 
axial heat conduction within the storage unit eliminated, 
so that the two systems are compared on a common basis. 

The Fortran code of the computer model and notes on its 
use are to be found in Appendix D. 
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Results of Computer-Model-Predicted TESS Designs 


The purpose of the dual heat exchanger TESS design exer- 
cise was to produce a cost optimized design of the dual 
heat exchanger TESS for comparison with the cost opti- 
mized design of the single heat exchanger TESS reported 
in NASA CR-159465. The single heat exchanger TESS re- 
quires the use of an external heat exchanger to vapor- 
ize and superheat the power fluid, whereas in the dual 
heat exchanger TESS this function is incorporated in- 
ternally and the external heat exchanger is eliminated. 
Therefore, the cost of the dual heat exchanger design 
is credited with the cost of the external heat exchanger. 
This is done in two steps j that is, the specific costs 
of the two designs are compared with, and without, the 
external heat exchanger credit. Additional savings in 
operating costs in the solar-electric system which accrue 
to the dual heat exchanger design are pointed out, but 
not credited. 

The specifications of the TESS design are as follows: 

Storage capacity is 3-1 x 10^ kJ when oper- 
ating over the range 516 K to 584 K. Heat 
is charged into storage medium (Thermkeep) 
by the non-phase change heat transfer fluid^ 

Caloria HT-43 aF a maximum rate of 1.8 x 10° 
kJ/h at an entering temperature of 584 K + 2 K. 

Heat is discharged from storage by means of 
toluene which enters as liquid at 473 K, and 
leaves as superheated vapor at 582 K + 2 K 
and 1,862 kPa. The maximum heat discharge 
rate is 1.0 x 10° kJ/h. The TESS includes 
an auxiliary heater which supplies auxiliary 
heat to raise the temperature of the toluene 
vapor to the specified temperature as re- 
quired. The Figure Of Merit (FOM) is the 
fraction of the heat delivered which comes 
from storage. 

Using the computer program described in Appendix D (Use of 
Computer Program) , trial designs were run on the Solar 
Cycle running through consecutive cycles to produce sta- 
bilization. These were chosen for convenience to have 
power fluid and heat transfer fluid heat exchangers of 
equal area, the inside heat exchanger surfaces of each 
varying from 37*2 m^ to 148.6 wr (400 to 1600 ft ; in 
increments of 18.6 m^ (200 ft^). Note that these are 
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the inside areas of each of the two heat exchangers, 
Thermkeep masses were varied in increments of 1 Mg 
between 10 and 25 Mg. A limitation imposed upon all de- 
signs was that the TESS always furnish power fluid of 
quality 1.0 vapor at the outlet; if it fell below this, 
the program was terminated. 

For a given heat exchanger design, a minimum Thermkeep 
mass was found below which the outlet vapor quality 
dropped below 1.0. For each design a storage specific 
cost and a Figure of Merit (FOM) were calculated. 

For comparison purposes the single-tube heat exchanger 
design program previously developed (NASA CR-159465) 
was rerun with several program changes to provide the 
same operating conditions as in the dual heat exchanger 
design; viz. axial conduction and insulation loss were 
eliminated, the Thermkeep properties were taken from 
temperature-enthalpy Table 10 (see Part II), and a re- 
turn temperature was calculated for the heat transfer 
fluid flow from the TESS combined with that from the 
external power fluid boiler, the results of which are 
shown on Figure 10. The Fortran code of the modified 
program is given in Appendix E. As a result of the modi- 
fications, the FOM rose from .88 to .92, and the specific 
cost decreased from $11.20 to $9.57/MJ, as compared with 
the values given in NASA CR-159465. 

The cost of each unit was computed using identical mater- 
ial and fabrication costs, tube wall thickness (.0004 m) , 
insulation thickness ( . 6 09 6 m) , aspect ratio (1.0), and 
vessel wall and bottom thickness (.0127 m and .0254 m, 
respectively) . 

Input data consist of the cost of material for each compo- 
nent and costs for fabrication on the basis of weight and 
surface area. Material costs were those at the time of 
construction of the experimental model in 1977. These 
fabrication costs are based upon the actual manufacturing 
cost of two experimental prototype "Thermbank" electric 
water heaters constructed at that time of the type de- 
scribed in the background section of the Introduction. 
Therefore, they do not include the cost benefits of the 
"learning curve" which results from converting prototype 
designs to manufacturing designs nor the cost reductions 
which result from quantity manufacturing. 
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Output costs incorporate fabrication factors based on 
weight, surface area, and the "six-tenths rule" whereby 
fabrication costs increase as the 0.6 power of the size 
ratio. These methods are generally accepted for making 
approximate estimates of costs of industrial equipment 
such as heat exchangers and vessels. 

Graphs were made of the parameters against specific cost 
Figure 11 shows the FOM against specific cost for several 
heat exchanger surface areas. The decreasing FOM seen 
in the 37-2 m 2 (400 ft 2 ) and 55-7 m 2 (600 ft 2 ) lines 
represent smaller Thermkeep masses, and these lines end 
when a smaller mass no longer produces quality 1.0 vapor 

Figure 12 shows the FOM against specific cost at con- 
stant Thermkeep masses. Here there is great improve- 
ment in FOM for relatively small cost increases. The re 
Iationships between the two graphs are plotted in Figure 
13 and l 1 * where the surface area and Thermkeep mass are 
varied at constant FOM. The Reference Design point for 
the single heat exchanger TESS is also shown on these 
figures . 

The choice of the Reference Design point is somewhat 
arbitrary, its choice depending upon the relative 
importance of specific cost and FOM. Figure 11 shows 
that an FOM as high as .99 can be achieved at a specific 
cost below that of the single heat exchanger TESS. 
Therefore — somewhat arbitrarily — a Reference Design 
which betters the single heat exchanger TESS in both 2 
specific cost and performance (FOM) is chosen at 55.7 m 
( 60Q ft 2 ) of each heat exchanger, 13 Mg of Thermkeep, 
and .95 FOM. It should be noted that, while the abso- 
lute specific costs are subject to the limitations 
noted above, the relative specific costs and perform- 
ances are realistic. The characteristics of the Ref- 
erence Designs are shown in Table V. 

The specific costs shown on the figures are computed on 
the same basis for the single, and dual, heat exchanger 
designs. However, as noted above, the dual heat ex- 
changer eliminates the need for an external toluene 
boiler. A credit for the boiler was estimated as fol- 
lows. A communication from the Sandia Laboratories at 
Albuquerque (Ref. 3) showed that the cost of the boiler 
in the solar total energy system which the TESS is de- 
signed to match, cost $17,000 in 1974-1975, and would 
cost approximately $25,000 to replace in 1980. 
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Storage specific cost, $/MJ 
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Figure 11. Specific cost vs. Figure of Merit for 
various heat exchanger areas. 
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Storage specific cost, $/MJ 
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Figure 12. Specific cost vs. Figure of Merit for 
various Thermkeep masses . 




TES specific cost, $/MJ 
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Figure 13. Specific cost vs. heat exchanger 
area for various Figures of Merit 


59 




60 


Figure lH. Specific cost vs. Thermkeep mass for 
various Figures of Merit. 


TABLE V 

CHARACTERISTICS OF REFERENCE DESIGN 


Cost Calculation Input Data 

Dual Heat 
Exchanger 
Design Point 

Single Heat 
Exchanger 
Comparison Point 

Shroud material cost, $/kg 

• 509 

• 509 

Shroud side fabrication cost, $/m^ 

10. l4 

10.14 

Shroud side fabrication cost, $/kg 

1.038 

1.038 

Shroud top fabrication cost, $/m^ 

77.35 

77.35 

Shroud top fabrication cost, $/kg 

7.92 

7.92 

Shroud bottom fabrication cost, $/m^ 

99.84 

99.84 

Shroud bottom fabrication cost, $/kg 

7.67 , 

7.67 

Shroud top thickness , m 

1.524x10 3 

1.524x10-3 

Shroud side thickness, m 

1.2192x10-3 

1.2192x10-3 

Shroud bottom thickness, m 

3.4036x10-3 

3.4036x10 -3 

Vessel wall material cost, $/kg 

.4718 

.4718 

Vessel end material cost, $/kg 

.4806 

.4806 

Vessel fabrication cost, $/m^ 

615.47 

615.47 

Vessel fabrication cost, $/kg 

2.447 

2.447 

Tube material cost coefficient 

.028 

.028 

Tube material cost, $/kg 

3.9091 

5.1998 

Tube fabrication cost, $/kg 

11.9 

11.9 

Tube fabrication cost, $/m2 

964.3398 

964.3398 

Storage material cost , $/kg 

.441 

.441 

Insulation cost, $/kg 

.2866 

.2866 


Cost Calculation Output Data 


Tube material cost , $ 

1479- 

2013. 

Total tube cost, $ 

14331. 

15190. 

Storage material cost , $ 

5733. 

7938. 

Insulation cost , $ 

992. 

1173. 

Shroud side material cost, $ 

173. 

199. 

Shroud side total cost, $ 

336. 

376. 

Shroud bottom material cost , $ 

120. 

138. 

Shroud bottom total cost , $ 

582. 

64o. 

Shroud top material cost, $ 

54. 

62. 

Shroud top total cost , $ 

315. 

345. 

Shroud total material cost , $ 

346. 

398. 

Shroud total cost , $ 

1232. 

1361. 

Vessel material cost, $ 

1321. 

1644. 

Vessel total cost , $ 

3271. 

3929. 

Total materials cost for TES, $ 

9872. 

13166. 

Total fabrication cost for TES, $ 

15687. 

16426. 

TES total cost, $ 

25559. 

29592. 

% OF TOTAL COST, HEAT EXCHANGER 
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This represents an annual cost increase over the six 
years of 6 . 6 %. On this hasis the cost in 1977, 
the year in which TESS comparative costs are based, 
would have heen $20,600. . This was subtracted from the 
computed TESS cost in order to obtain a Storage Spe- 
cific Cost After Boiler Credit. 

An additional benefit would accrue to the dual heat ex- 
changer TESS in the Sandia solar total energy system 
resulting from the elimination toluene holler; namely, 
the heat loss from the boiler during the nighttime sys- 
tem shutdown, especially in winter (Ref. -4). Since this 
represents an operating cost benefit, no credit is 
taken in the capital costs reported here. 

The side-by-side tube concept was chosen for this study 
because it is simple to represent mathematically, and a 
method of constructing such tubing can be easily visua- 
lized. For small diameters (.004 - .03 m) , the cost per 
unit weight of tubing decreases as diameter increases 
and this is reflected in the heat exchanger costs of the 
Reference Design. A double tube consisting of two identi- 
cal tubes was considered to have a material cost just 
twice that of a single tube. These factors may unfairly 
favor the dual tube costs since they are based upon a 
tubing configuration which is not known to be readily 
available. However, it seems probable that such tubing 
can readily be manufactured. For example, two tubes 
might be passed side-by-side through a continuous MIG 
(metal-inert gas) welding machine which lays a fillet 
between them, or an electric current could be applied 
across them to cause a continuous line weld between them. 

Another concept which might apply would use two coaxial 
tubes with an annulus between, plus a single tube. One 
fluid would flow through the annulus, while the other 
fluid would flow through the inner coaxial tube and 
through the single tube. This system provides for direct 
heat exchange between the two fluids, and between each 
fluid and the storage medium. Coaxial tube configura- 
tions are presently available; designs exist with fins 
and wire in the annulus to promote turbulence or direct 
heat transfer between the inner and outer walls. 

The computation of the single and dual heat exchanger 
costs on the same material and fabrication cost basis 
may underestimate the cost of the dual heat exchanger. 

In order to show the latitude for cost increase, the 
cost of the dual heat exchanger was modified in the 
following three cases: 
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Case 1 


The material cost of tuhing of the dual heat 
exchanger was fixed at twice that of the 
single tube exchanger f$10.IJ/kg vs. $5. 2/kg) 
to eliminate the cost advantage of larger 
diameters. Then the fabrication costs were 
allowed to increase 22 % until the storage 
specific cost of the dual heat exchanger TESS 
was the same as that of the single tube heat 
exchanger TESS, viz. $9-.57/mJ. 

Case 2 


The material cost was kept the same ($3- 9/kg), 
and it was found that the fabrication costs 
could be allowed to increase bl% before the 
storage specific cost of the single tube heat 
exchanger TESS was exceeded. 

Case 3 


The material cost was kept the same but the 
fabrication costs were doubled. This pro- 
duced a storage specific cost of $11.91/mJ for 
the dual heat exchanger TESS. 

The results of these comparisons are shown in Tahle VI. 

Table VII shows the percentage of total cost of the three 
cases of the dual heat exchanger TESS described above, and 
for the single heat exchanger comparison case. 

The characteristics of the Reference Designs for the dual 
and single heat exchanger are summarized on Table VIII 
where the Case 3 heat exchanger cost has been assumed for 
the dual heat exchanger TESS in order to present the re- 
sults on what is believed to be a conservative basis. On 
this basis the dual heat exchanger TESS has a specific 
cost which is higher than the single heat exchanger TESS 
by 25 % before the boiler credit, but which is only 58% 
of the single heat exchanger TESS cost after boiler credit. 
This cost saving is achieved while the performance as 
measured by the Figure of Merit is improved from .92 to 
• 95. 
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TABLE VI 


COMPARISON OF THREE METHODS OF CALCULATING 
HEAT EXCHANGER COST FOR DUAL HEAT EXCHANGER TESS 


Case 1 Case 2 Case 3 

Cost Calculation Input Data: 


Tube Material Cost Coefficient 

.0745 

.028 

.028 

Tube Material Cost, $/kg 

10.401 

3.909 

3.909 

Tube Fabrication Cost, $/kg 

14 . 518 

16.815 

23.8 

Tube Fabrication Cost $/m^ 

1176 

1362 

1929 


Cost Calciilation Output: 


Tube Material Cost, $ 

3934 

1479 

1479 

Total Tube Cost , $ 

19614 

19639 

27183 

Total TESS Materials Cost $ 

12327 

9872 

9872 

Total TESS Fabrication Cost $ 

18515 

20995 

28540 

Total TESS Cost $ 

30842 

30867 

38411 

Storage Specific Cost, $/MJ: 

Before boiler credit 

9-57 

9-57 

11.91 

After boiler credit 

3.18 

3.18 

5.52 



TABLE YI1 


PER CENT OF TOTAL COST OF TESS 



Dual Heat Exchanger 
Case 1 and Case 2 Case 3 

Single 

Heat 

Exchanger 

Thermkeep 

19 

15 

27 

Heat exchanger 

6k 

71 

51 

Containment 

11 

9 

13 

Thermal insulation 

3 

3 

k 

Outer shroud 

k 

3 

5 
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TABLE VIII 


SUMMARY OF DUAL AND SINGLE HEAT EXCHANGER 
TESS DESIGN POINTS 



Dual 

Single 

Tank height and diameter, m 

2.11 

2.35 

System height and diameter, m 

3.4 

3.6 

Quantity of Thermkeep, kg 

13,000 

18,000 

Heat exchanger: 



Tube inside diameter, cm 
Number of tubes 
Total length, m 

.716 

150 (prs ) 
5588 

0.457 

280 

7761 

Figure of merit 

• 95 

.92 

Total cost $ 

38411 

29592 

Specific cost, $/MJ 

11.91 

9-57 

Boiler credit $ 

20600 

None 

Net specific cost ? $/MJ 

5.52 

9-57 


66 



CONCLUSIONS 


A program of experimentation and analysis has been 
carried out to continue and extend the work reported 
in NASA CR-159465, on the development of a phas ex- 
change thermal energy storage system (TESS] using modi- 
fied anhydrous sodium hydroxide for solar electric 
power generation. 

The phase-change medium used comprises commercial grade 
anhydrous NaOH modified by the addition of 8 % (wt) of 
NaNO^ and 0.2$ MnOg* 

The TESS studied in the prior work utilizes a single 
heat exchanger and a non-phase-change heat transfer 
fluid, and is designed to meet the following operating 
conditions : 

Storage capacity is 3-1 x 10^ kJ operating 
over the range 516 K to 584 K. Heat is 
charged into the storage medium at a maxi- 
mum rate of 1.8 x 10° kj/hr, at a tempera- 
ture of 584 K + 2 K. Heat is discharged g 
from storage at a maximum rate of 1.0 x 10° 
kJ/hr at a temperature of 582 K + 2 K. 

The present program consisted of three parts: 

Part I extends the prior work by means of 
chemical analyses of the storage medium, 
and examination of the heat exchanger of 
the one-tenth scale experimental model 
constructed and tested in the prior program. 

Part II utilized the physical model of the 
prior program to obtain experimental data 
with a different heat transfer fluid, 

Caloria HT-43, for further verification of 
a computer model of the prior program which 
used Therminol-66 . 

Part III included the development of a new 
computer model which introduced a second 
heat exchanger into the TESS for evapora- 
ting and superheating the power fluid for a 
Rankine cycle power generator. This was 
used to produce a cost optimize reference 
design of a TESS for comparison with the 
single heat exchanger design of the prior 
work. 
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No chemical deterioration of the medium was found. Some 
physical segregation of components was observed, which 
is attributed to the fact that the medium is a non- 
eutectic mixture, but no corresponding deterioration in 
performance was noted. Physical segregation can be 
reversed by periodic complete melting of the medium. 

No significant physical deterioration of the heat ex- 
changer was observed. 

The dual heat exchanger TESS was found to be superior to 
the single heat exchanger design in performance and in 
capital cost, due to savings resulting from the elimi- 
nation of an external heat exchanger. Additional ad- 
vantages result from system operating cost savings. 

Auxiliary heating is used to raise the temperature of 
the power fluid vapor to the specified temperature, and 
a figure of merit is defined as the ratio of the heat 
from storage to the sum of the heat from storage plus the 
auxiliary heat. This is used as the measure of performance. 

The design comparison is summarized as follows: 


Dual Heat Single Heat 
Exchangers Exchangers 


Quantity of storage medium, kg 
Heat exchanger tubing 

13,000 

Inside diameter, cm 

0.716 

Total length, m 

5588 

Total cost 

$38,411 

Specific cost, $/MJ 

11.91 

Boiler credit 

$ 20,600 

Net specific cost, $/MJ 

5.52 


18,000 

0.457 

7761 

$29,592 

9.57 

- 0 - 

9.57 


Figure of Merit 


•95 .92 
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■TECHNICAL REPORT- 
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CLIENT : Comstock 6 Nescott, Inc. 

CASE NO : 11833 

3. Colorimetric Method, using spectrophotometer, Beckman DB 
for sodium nitrite. 

4. Official Methods of Analysis of the Association of Official 
Analytical Chemists. 

5. Scott's Standard Methods of Chemical Analysis. 


RESULTS OF ANALYSIS: 




Sodium 

Nitrate 

Sodium 

Nitrite 

Sodium 

Carbonate 

Sodium 

Chloride 

Manganese 

Dioxide 

Sample 

#1 

8.467 

0.0306 

Percent 

2.62 

2.03 

0.172 

Sample 

#2 

7.477 

0.0307 

1.45 

1.98 

0.171 

Sample 

#3 

4.652 

0.0020 

2.85 

1.73 

0.101 

Sample 

#4 

5.278 

0.0021 

1.48 

1.35 

0.085 

Sample 

#5 

5.041 

0.0023 

1.55 

1.47 

0.113 

Sample 

#6 

5.326 

0.0022 

1.43 

1.63 

0.110 

Sample 

#7 

4.503 

0.0019 

1.38 

1.68 

0.083 

Sample 

#8 

4.377 

0.0023 

1.35 

1.85 

0.110 

Sample 

#9 

5.623 

0.0023 

2.89 

1.65 

0.102 

Sample 

#10 

4.878 

0.0026 

1.36 

1.76 

0.095 

Sample 

#11 

2.791 

0.0017 

1.44 

1.53 

0.042 

Sample 

#13 

2.707 

0.0021 

1.40 

1.76 

0.057 

Sample 

#14 

2.546 

0.0032 

1.32 

1.25 

0.053 

Sample 

#15 

5.772 

0.0098 

2.88 

1.58 

0.137 

Sample 

#16 

5.138 

0.0095 

2.90 

1.43 

0.172 
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Chemist 
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APPENDIX B 


HYSTERESIS STUDY GRAPHS 
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SYMBOLS 


Thermkeep temperatures: 


o = Shallow well temperatures vs. 
location In tank 


A = Deep well temperatures vs , 
location In tank 


+ = Computer model predictions 
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Heat transfer fluid (HT-43 and T-66) temperatures: 

o = Experimental data 

= Computer model predictions 

Thermkeep temperatures: 

o = Shallow well temperature vs. location 
in tank 

A = Deep well temperature vs. locations in 
tank 

□ = Tank surface temperature vs. location 

in tank 

+ = Computer model prediction 
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Figure C-27 • Thermkeep temperature profile, charge test 215-2 (l.l gpm) , after 6 hours 
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Figure C-28, Thermkeep temperature profile, cyclic test 218 Cl,Q gpm)_, starting profile. 
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Figure C-29. Thermkeep temperature profile, cyclic test 2l8 (1.9 gpm , at end of first discharge phase. 
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Figure C-31. Thermkeep temperature profile, cyclic test 218 ( 1 . 9 gram), at end of second discharge phase 
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Figure C-32, Theimkeep temperature profile, cyclic test 2l8 (.1.9 gpml, at end of second charge phase. 
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Figure C-33. Thermkeep temperature profile, cyclic test 2l8 (1.9 gpm), at end of third discharge phase. 
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Figure C-3^. Thermkeep temperature profile, cyclic test 218 Cl. 2. gpml, at end of third charge phase 
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Figure C-35- Thermkeep temperature profile, cyclic test 218 (1.9 gpm) , at end of fourth discharge phase 
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Figure C-38. Thermkeep temperature profile, discharge test 212_ (J .7 gpml, after 2 hours 
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Figure C-39* Thermkeep temperature profile, discharge test 219 (7.7 gpm) , after H hours. 
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Figure C-4l. Thermkeep temperature profile, solar cycle test 223, at end of Phase 1 discharge, cycle 1. 
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Figure C-43. Thermkeep temperature profile, solar cycle test 223, at end of phase 2 discharge, Cycle 1 
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Figure C-L4. Therakeep temperature profile, solar cycle test 223, starting profile, cycle 2. 
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Figure C-H 5 - Thermkeep temperature profile, solar cycle test 223 , at end of phase 1 discharge, cycle 2 
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Figure C— 46, Thermkeep temperature profile, solar cycle test 223 ? at end of cycle 2 charge 
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Figure C-UT. Thermkeep temperature profile, solar cycle test 223, at end of phase 2 discharge, cycle 2 
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Figure C-U8. Thermkeep temperature profile, solar cycle test 223, at end of overnight idle after cycle 2 
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Figure C-51. T-66 outlet temperature, charge test 004 (3.4 gpm) 
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Figure C-52. T -66 outlet temperature, discharge test OlU (JL,Q gpml. 
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Figure C-53. T -66 outlet temperature, charge test 015 (l.O gpm). 
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Figure C-5U. T -66 outlet temperature, cyclic test 0l8 ( 2.0 gpm) , first cycle. 
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Figure C-60. Thermkeep temperature profile, cyclic test 018 (_2.0 gpm), starting profile. 
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Figure C-6l. Thermkeep temperature profile, cyclic test 018 (2.0 gprn), at end of first discharge phase. 
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Figure C-62, Thennkeep temperature profile, cyclic test Ol8 C2,Q gpmi, at end of first charge phase 
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Figure C-63. Thermkeep temperature profile, cyclic test 018 (2.0 gprn) , at end of second discharge phase. 
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Figure C-6U, Thermkeep temperature profile, cyclic test Ql8 C2,Q gpml, at end of second c har ge phase 
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Figure C- 65 . Thermkeep temperature profile, cyclic test 018 (2.0) gpm) , at end of third discharge phase. 
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Figure C-66. Thermkeep temperature profile, cyclic test Ol8 (2.0 gpm), at end of third charge phase 
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Figure C-68. Thermkeep temperature profile, solar cycle test 023, at end of phase 1 discharge, cycle 1. 
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Figure C-69. Thermkeep temperature profile, solar cycle test 023, at end of cycle 1 charge 
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Figure C-70. Thermkeep temperature profile, solar cycle test 023, at end of phase 2 discharge, cycle 1. 
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Figure C-71. Thermkeep temperature profile, solar cycle test 023, starting profile, cycle 2 
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Figure C-72. Themkeep temperature profile, solar cycle test 023, at end of phase 1 discharge, cycle 2. 
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Figure C-75. Thermkeep temperature profile, solar cycle test 023, at 32 hours 




APPENDIX D 


USE OF COMPUTER PROGRAM 


Input Data 


In order to run a test case with the FORTRAN program, 
data must be entered into the DATA section in the proper 
format. Each piece of data which is read by the program 
is then identified and printed in the output section of 
every run to provide easy checking of all parameters. 
There are three main groups of data: 


1. Main Program Data 


A table of solar collector output vs. time, 
fluid characteristics, and program control data 
(Fig. D-l ) . These are self-explanatory except 
for the following: 

ALLOWABLE POWER FLUID OUTLET TEMP is an arbi- 
trary temperature, lower than the desired out- 
let temperature, used to calculate the amount 
of auxiliary heating necessary if the design 
criteria were relaxed to this point. 

STABILIZATION PERCENTAGE is the percentage dif- 
ference between storage input and output over 
a complete cycle which must be achieved before 
the program will stop. At this point the unit 
is assumed to be repeating itself each day. 

NUMBER OF DAILY CYCLES is the maximum number of 
cycles the program will be allowed to run. 


2. Element Analysis Subroutine Data (Fig. D-l) 


Definitions : 

p 

AEL: The total area in m of inside surface 

of the transfer fluid heat exchanger 
in contact with storage medium. 
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2 

AEP: The total area In m of Inside surface 

of the power fluid heat exchanger in 
contact with storage medium. 

2 

ADL: The total area in m of common direct 

contact between the power fluid and 
transfer fluid heat exchangers. 

DEPF: The effective diameter in meters of 

one double tube of heat exchanger 
upon which storage medium solidifies 
(the diameter of a circle whose 
perimeter is that of a cross section 
of a touching power fluid and transfer 
fluid tube ) . 

XLEFF: The effective length in meters of one 

tube of heat exchanger. 

XNT: The number of heat exchanger tubes. 

RHOM: The density in kg/m J of storage material. 

HTCS: The heat transfer coefficient of suh- 

cooled power fluid in kW/m 2 -°K. 

HTCV: The heat transfer coefficient of power 

fluid vapor in kW/m 2 -°K. 


Thermkeep Properties Subroutine Data 


D-l ) 


In addition to the DATA section, two functions are 
modified directly for each run: HL(wL) and HT(X). 

HL(wL) computes the heat transfer coefficient of 
the transfer fluid as a function of flow rate in 

kJ/sec-in -°K . For this function the tube size, the 
fluid properties and a helix diameter must be known. 
HT(X) computes the heat transfer coefficient of 
vaporizing power fluid as a function of quality in 

kJ/sec-m 2 -°K. 
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DUAL HEAT EXCHANGER - HEAT OF FUSION THTRMAL STORAGE ANALYSIS 
TADLF OF SOLAR COLLECTOR OUTPUT VS. TIME 
RATE IN KW TIME IN SECONDS 


c.c 

O.C 

c.c 

8800.0 

2 5 C . C f 0 C 

13900.0 

323. COOO 

15600.0 

358. PCOO 

16600. C 

432.0000 

18900. C 

4 4 e . C C 0 0 

19900 .C 

471.0000 

21600. C 

480.0C00 

22600.0 

498.CC GO 

25100.0 

f CC. OOOC 

259 00. 0 

493.0000 

27800.0 

470,0 COO 

30300.0 

442. COCO 

32700.0 

3 1 2. 0000 

35600.0 

250. COOO 

39 0 OC . C 

r . o 

43200.0 


MAIN PROGRAM DATA 


TRANSFER FLUID SPECIFIC Hr A T . KJ/KG/K 
TRANSFER FLUID INLET Tf fPE R A T URL . K 
PPHTR FLUID SP .HT .-LI QUID. KJ/KG/K 
POWER FLUID HFAT/VAP, K J/KC 
POWER FLUID SP.HT.- VAPOP. KJ/KG/K 
P 0 W T R FLUID SATURATION TLMPLRATURE, K 
POWER FLUID INLET TEMPERATURE. K 
DESIRED PDwtR FLUID DUTU'T TEMP., K 
AL LOW AL j LE FOWtR FLUID EjUTLET TEMP. » N 
POWER FLUID FLOW RATE, KG/SEC 
INSOLATION STARTING T 1 MT » SEC 
BEGIN NET CH A PGE . SEC 
END NTT CHARGE, SEC 
THFRMKEEP SPFC1EIC HEAT, KJ/KG/K 
NUMBER OF ELEMLNTS FOR ANALYSIS 
AMOUNT OF STORAGE MATERIAL. KG 
MINIMUM TIME INCREMENT. SFC 
MAXIMUM TIME INCREMENT, SEC 
MAXIMUM TRANSFER FLUID FLOW, KG/SIC 
STAT 1 L I7AT ION PE RCFNT ACT 
NUMBER OF DAILY CYCLES 


r. 5?sc 

5 8 4 .00 
2.260000 
217 .03899 fi 
1 .957000 

530.00 
A 7 3 . 0 0 
58A .00 
5B0.0C 

0.55342 

8500.00 

1 3500.00 
39000. CO 

2.5000 

1C. 

10000.00 
15.00 

100.00 
5. COOO 

3 .00 
6 


ACP = 
PE F F = 
RHOM = 


A 9 . 5 4 3 3 AFL= 49.5483 APL= 6.1935 

0.013009 XI F FF = 18.628 XNT= 150.000 

1802.080 HTCS= f.39674C HTCV= 0.976010 


TABLE OF TE 

MPER ATURF ,CN TH4LPY, 

AND SOLID FRACTION 

T ( I ) 

HE 1 ) 

Z ( 1 ) 

372. CO 

- 1 3 .* COOC 

1 .000000 

505.60 

259.64990 

1.000000 

528.00 

337.36011 

1 .000000 

53C. 00 

344 . 30005 

0.760000 

535.00 

-»61 .64990 

0.736200 

541.63 

4 09.8 60 1 1 

0.670000 

544.89 

4 33 .6 0 0 1 0 

0.615000 

549.17 

464 .69995 

C .59CU00 

557.86 

628.0 COOC 

C .4 OOC 00 

560.00 

5 4 3.5 70 0 7 

0.368700 

563. OC 

600 .25000 

0.347500 

573.00 

672 . f 00 0 0 

0.2 9 5' 00 0 

583. OC 

7 11 . CCOOO 

0 .0 

643.00 

8 36 .8 50 1 0 

C.O 


OF TMET M KF E P 


Figure D-l. Dual heat exchanger program input data. 
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Input Data Selection 


Much of the input data does not change from run to run. 
The insolation table, Thermkeep properties, and all the 
main program data remain constant except for the Number 
of Elements for analysis (normally 10 for general use, 
and 25 for more accurate final designs) and the amount 
of storage material. The insolation table is a copy of 
the "Albuquerque Winter Mean Solar Day" modified to 12 
hours by omitting a turbine overnight idle period. The 
Thermkeep properties are "Table 10" for Thermkeep de- 
termined by the work described in Part II. The remain- 
ing data applies to toluene as the power fluid and the 
system characteristics as required by this contract. 

The maximum and minimum time increments were determined 
by trial and error during the early stages of program 
development . 

The element analysis data were chosen with the help of 
a desk top computer. First, desired heat exchanger 
sizes, areas, and tube diameters are arbitrarily selec- 
ted. Then the number and length of the tubes are varied 
until a combination is found which produces about a 10 
psi drop in the transfer fluid exchanger and a reason- 
able drop ( 20-30 psi) in the power fluid exchanger. 

Then an amount of storage material is arbitrarily se- 
lected. If the height-to-diameter ratio is specified 
(normally 1.0), this produces a helix diameter of the 
tubing. From the known mass flow of power fluid and 
heat exchanger geometry, the heat transfer coefficients 
of the power fluid are now calculated. 

The effective diameter is calculated by using the diame- 
ter of a tube whose perimeter is the same as the sum of 
the two fluid tubes’ perimeters minus the section which 
they have in common. This report uses all equally-sized 
tubes for both power and transfer fluids and the common 
area was 11$, or 40° of arc on each tube in direct con- 
tact (Figure D-2). A small addition to the perimeter was 
made for weld fillet between the tubes. 

Function HL(wL) receives the transfer fluid flow rate 
from other parts of the program and returns a heat trans- 
fer coefficient. This function is modified for each run 
because it uses the tube inside diameter and helix 
diameter which change with each case. The helix diameter 
is calculated as the diameter of one heat exchanger coil 
such that one-half of the storage material is enclosed 
by the coils and one-half lies outside them. 




From the heat exchanger geometry and transfer fluid prop- 
erties, the following coefficients are determined: 


CRE: 


RECRIT: 


CHTC : 


QX: 


HTCI : 


A multiplier of the mass flow rate 
which will give the Reynolds number 

CRE = d^/C fluid area x fluid viscosity) 

where d^ is the tube inside diameter. 

The critical Reynolds number 

RECRIT = 2100 Hl + 12 (d + d h ) 0,5 

where d, is the helix diameter, 
h 

A heat transfer coefficient factor: 

CHTC = (.022 x C f x Pr* i| )/d t 


where C„ is the fluid thermal conduc- 
tivity . 


A dimensionless factor used In the cal- 
culation of the effect of the helix 
diameter on heat transfer in laminar 
flow : 



The heat transfer coefficient for turbu- 
lent flow: 


HTCI 


CHTC x Re 


8 
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Function HT(X) receives the power fluid quality and re- 
turns a heat transfer coefficient. The relationship 
for toluene is roughly a sigmoidal shaped curve, of 
which a linear approximation is a good average. There- 
fore, during this contract we used 


HT( X) = X x (HTCV - HTCS ) + HTCS 


Output Data 


For each case run, a back-up calculation print-out was 
made. An example is shown in Figure D-3- Cost calcula- 
tions are based on the actual cost of a prototype stor- 
age unit, and integrate manufacturing and material costs 
both on a surface area and weight basis. The only dif- 
ference from case to case is in the tube material cost, 
where there is a lower cost by weight as tubing gets 
larger. 

A dollar amount for each category is then printed out, 
along with a unit total cost. As a point of interest, 
the percentage of total cost of the heat exchanger, insu- 
lation, cover, tank and storage material are also listed. 

Next, the geometry and weight input data are tabulated, 
and the resulting dimensions and weights are listed 
(Figure D-4 ) . 

Finally, the heat exchanger dimensions are calculated 
and listed along with Reynolds numbers, heat transfer 
coefficients, and pressure drops. These outputs allow 
us to check total surface area and to vary the number, 
length, and diameter of tubes in order to produce desired 
pressure drops. With toluene, we always found that it 
was in laminar flow in the liquid region, turbulent flow 
in the vapor region, and the vast majority of pressure 
drop was due to the vaporization acceleration. Also 
calculated were the power fluid pressure drop assuming 
that it was liquid through the entire tube, and vapor 
pressure drop assuming it was vapor for 90 % of the 
tube — both of these cases were assumed extremes. 



COST CAI.CUI.Af JON INPUT HAT A 


3IIR0UH MATERIAL COS I , i/KU. .009 

SI IROlJTi S I I>L FABRICATION COST, i/SO.M, 1 (1 , 1 1| 

S|lROUl.i-S.Ll.i£.. I A.UJUCAT J U N ililil _i > LUu J . 0 M 1 

si ikoun I or 1 r auk i ca t ton cos i , i/so , m . it. 35 

HUROUII TOR FABRICATION COST, i/KG, t .9? 

SH..RQUU _BO J . LQM-.F.ADR 1 GAUjQiL -CQliU-J'/SQ . Ji , . 99,04 

SHROUD BOTTOM FABRICATION COST , T/KG . t . 6 t 

SI IROllli T OP fl 1 1 CR NT SS , M , I . 5?4 I- - 3 

iiLLROilli SI Hi 1 II I SS , H . 1 . ? I RIM y , 

SI IROUB BOTTOM Till CKNI SS, M. 3 . 4036F- 3 

VrSSIl. WAI.L MATLRJAI, COST. i/ICO. ,4710 

VI--SSI-I. I . N I.i MATT RIAL. COST, t/ICO .4006 

VI SSI I , T ABRTLA T LON ('OS I , 1 ./CO . M , 610,4 t 

VKSSI.I , I ABR I f'A T I ON COS T , i/ICO . ? . 4 4 7 

TORI MATERIAL COST C OI FI 1C J | N"l _ . 0 P t J 

7 1 1 BL- MAT L R I Al . COST, i /K 0 , 6 , 9 0 V 1 

TI.JBL I ABRJ CAT I ON COST . i/l<0. 11.9 

Tl I Til" FAB R T CA T I ON COS T . i/S Q . M ■ . 964 , 3398 

STORAOT MATIRIAC COS i , i/|<0. ,44 1 

JNSUi.AT T(TN COST . i/ICO. .9066 


COS T CALCULA7 l ON 00 I PU I HA T A 

film MATLR.I.AI ■ COST, i 

TOTAL TUBb COST, i 
S FORAGE MA Tl- R I Al . COS T , i 

iNSULAUON-CUij U ± 

SI IROl Hi 3 J HI MA I L R I AL, COST , i 
SIIROOH SJUI TOTAL COOT , T> 

SHROUD BQ I | QH MAT LR1 AL CO ST, _t 

SHROUD BOTTOM TOTAL COST, i 
SIIROOH TOP MATT- RI AL COST, t 

SIIROOH TOP TOTAL COST , 1> 

SIIROOH TOTAL MAT I RT Al, COS l , t 
SIIROOH TOTAL COST, i 

VI SSr-l, MAT I RIAL COST , j> 

VIISSLL TOTAL COST, 1> 

TOTAL MATER TALS COST FOR TIS, $ 

TOTAL FABRICATION COST T OR TIS, 

TFS TOTAL COST, T> 

% Ml- IHTAL COS) 


HEAT LXCIIANOI R 

6 1 

INSULA l ION 

4 

SHROUH 

0 

VL SSI; L 

1 ? 

TUI RMKFFR 

J 9 


Figure D-3. Cost calculations print-out example. 


1 4 79 , 

I 4 <3.1 . 
4 4 I 0 , 
073 , 

I 06 , 

6 0S , 

JOS, 

041 , 

4 S' . 

:> 93 . 

3 1 ? . 
114?, 

1113. 

?030, 

0 I 06 , 

1 04 00, 

?309 J , 
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pLOHE_ I RY AND WEIGHT INPUT BAJA 


TRANSFER - Ti l B I — 0 . H ■ y M ; 7 . 1 62 8 L -- 3 

T RANGIER TUHr I . D . , M. 6 , 30E -3 

I'ilUl R TIJUT O.H. , M. /. lAl'Hr 3 

POUI I -R TUB E: — L , H ■ , M . 6-r-60F-™3 

I DDF LENGTH, M. 10-620 

NO . OF T UHLS tOO 

I msi H. . AW I N TH rOK - N T 00 , - ti n ^64194 

Vh OOF L UAI .1. THICK Nl OS , M . .0127 

VFSOI I, r Nil THICKNESS, M. .0201 

T I I f M44! HI NS1 T - Y , KG/CU - .PL 4434)364143 

STEfc I. DLNOl TY, KO/CU.M. 7807 

FLUID DENSITY, KG/CU.M. 024 

INSULATION DENSITY, KO/OII.H. 140-~X2 

ASRLCI RATIO, HGf/DIA. 1 


AM ) . 01" 0(0 RAUF MAT F R I AL , KG. J 0 0 0 0 
GFOMF-FRY-- AND UE J Gl IT OUT PU T DAXA 


SHROUD D I AMFTCR, M. . 0.1 760 

VF0SFL - 4M A MF1T.R, M . L^lSUl 

TUBE HELIX D T AML" T I-. R , M . .111 S 

VI SSEI. HEIGHT , M. 1 .9310 

VFSSFL UF-l-nur-,— 4411, 23 3 A ■ 6 b 0 0 

ULIGHT or TUBES, ICO. 3 70.27 7 t 

INOUI.ATTUN UF I0IIT . KG. 3044 . 421 7 

OYOI1- M I I F I O UT, M . 4, -24118 

SYSTEM UFTGHf, KG. 16371. 4102 


POUt R TUBE 

ID “ 


.20 

IIFL1X D 1 AMT 

rr r " 


M.391 

NUMDI R 01 TUBI S 


ISO 

TUBI- l.l NOTH 

, 1 1 


61 . 1 1 US 

TUBE TNSinr 

ARL A, 

f 12 

60 ft . 

rirnr rNoiru 

ART" A. 

m2 

06 , 

L.TQUID 16 



7' 0 0 . 

VAPOR Re 



56S9II . 

CRTTTOAI. Ro 



01 1 3 . 

HTCL Btu/hr 

( (2 o 

F 

69 . 090 79 

H 1 Cl . k W/m2 

°l< 


. 396 74 

IITOV IMu/ln 

112 » 

F 

1 7 J . 94 907 

mrv uj/fl,*' 

4< 


. 976ft 1 

AP AruUp) <il ion 

PS 1 

10.3 

AP l i nil i (1 wh ole l ub e PSJ 

. 4 0 


AP v<ipo> ,90% lubf PS I I .*>/, 


AP T 66 whole lube POI 9.62 


Figure D-4 . Output data print-out example 



The FORTRAN program prints out all the input data as 
previously mentioned. Then the unit starting condi- 
tions are calculated and listed — the starting tempera 
ture, starting solid fraction, and total unit energy. 
The charging fluid flow rate and both fluid outlet 
temperatures are printed each half-hour of elapsed time 

A map of the unit showing the temperatures of storage 
medium and both fluids, along with fluid quality and 
solid fraction of the medium in each element or "slab" 
is printed out at two points in each daily cycle — at 
the beginning of net charge when the unit is at its 
lowest state of charge, and at the end of charge when 
the unit is at its highest state of charge. Also the 
total unit energy at these times is printed out, and 
the difference between them, once the unit is stabi- 
lized, is equal to the amount of energy actually output 
from storage, and does not include the heat transferred 
directly from the transfer to the power fluid. 

At the end of a daily cycle, a summary of the energies 
transferred is printed: 

The LOST COLLECTOR HEAT is the amount of 
heat available but unused because of maxi- 
mum flow conditions . 

The STORAGE OUTPUT is the total amount of 
energy which came out of the unit in the 
power fluid stream. 

The STORAGE INPUT is the amount of energy 
which entered the unit in the transfer 
fluid stream. 

The AUXILIARY INPUT is the amount of heat 
required to raise the power fluid outlet 
temperature to the Desired Outlet Tempera- 
ture . 

The DEFICIT REQUIREMENT is the amount of 
heat required to raise the power fluid out- 
let temperature to the Allowable Outlet 
Temperature . 

The EXCESS HEATING is the amount of heat 
which exited the unit above the desired 
outlet temperature. 
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The BOILER REQUIREMENT is the total daily 
boiler consumption. 

Also listed are the total collector output, the elapsed 
time, the number of time steps used during this cycle, 
and the cycle ending total energy. 

At this point, a new cycle is begun; or, if the unit 
has stabilized, the program stops. 
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COMPUTER CODE 


) . T I *.L ( 25 ) 


u u 

1C 1 

l'.r 

ICO 

10 4 

1C' 
K t 

rr i 
2C? 
Z 0 3 
7 C * 
< SC 
zot 
Z 0 7 
2 C o 
2CE. 
2 1 0 
2 1 1 
33 0 
33 1 
; i? 
2 1 3 
214 

1 1 f 

2 1 9 
3C 1 
1C c 

31 4 
317 

n c 

21 7 

? 1 P 
0 

32 1 
323 
32 4 
32 5 
32 C 
32 7 
32 P 
329 

Tir 

32 f- 
339 
24 (. 


1 ) 


f 1 


TAIN IKOf.UM Trrowj 

C OMM CK TP CUT .TLOUT , T 1 1 t . TL IN' . TSAT »CT *CL. CS.CV .CT • XNt .A ST OR , XOUT *TA 
>'R R R . T f .1 (i .NC 
r IF L Ni 1 ON Ql <2 

r r p p a t ( 1 2 ) 

F TP PAT ( F 1 C .4 ,r 1C 
FIPPAT(P1C.4 * 3F 1 
FCRVAT(4F10.2) 
rcpvAT<F ic.t »rnc.r,F 

F CRM/ T (2F 1 0. 2 ,F 1 C .4 ) 
f Cf MAT( 2F 1 0. 2 ) 
f rir wA T( 2f,x 

F C R m A T ( bX , 

F TP MAT < 5X . 

FT Ml A 1 ( NX , 


1 C . 4 . F 1 0 . C . 1 3 . F 1 0 . 2 ) 


F 1 G.4 «9X. T 1 3. 1 ) 

TRANSM R F Lu 1 F' SPECIFIC HEAT. FJ/KC/K 
PC Vi t r r L U H CP .HT. -LIQUID. KJ/KG/K 
P C *• t F FLUID HC A T /V A P » K.J/KG 
FTwri, FLl'ir S P . HT . - VAPOR , XJ/KG/K 
POET' FLUID SATURATION TEMPERATURE. K 
TPANSFC-? FLUID INLET TfMPCRATURC. ¥ 

[ C5 I R F p POvifF. FLU1C OUTLET TEMp., ¥ 
Pfv.fi, FLCU’ 1NCFT TEMPERATURE. K 
PC W f_ R FLUP c LO^ RATE. KG/SEC 
IN'SCLATIlN STARTING TOE. SEC 
ELGIN NET CHARGE. SEC 
END N F T CHARGf . S£ C 
AMOUNT OF STORAGE MATERIAL. KG 
THFRMKfrr EPFCIFIC HEAT. KJ/KG/K 
NUMPEP OF r Lf MCFiTS FOR ANALYSIS 
NL'HEFR OF DAILY CYCLES 

STAFIl izatiln pfr cent age 
• TArir OF SOLAR COLLECTOR' out rut 


rn m at ( sx . 
r CROAT ( ox . 

F CPM AT ( NX , 

ftrv at< bX , 

F or vat ( t X . 

F c R M A 7 ( r- X . 

F CRM A T ( OX , 

FCFMAT(5X, 

F LR M A T ( 0 X . 

FORMATE sx. 

F LR MA T { OX , 

FOR Mill bX, 

FORMA 7 (OX , 

F Lf VAT ( OX . 

FORMAT ( 1 9X 
r DRV AT ( 2 4 X 
F C PVATf 1 HI 
f PR.VATC 1H1 

set analysis 

F CRM AT < OX . 

F LRf'ATt OX . 

F C R V A T ( 0 X , 
f CRM A T ( /// 

F CRV AT < 1 OX 
F C R M A T < 1CX 
F CRV A T { 1CX 
F CR'M A T ( 10X 
FORMAT t 1 OX 
rcPMAT( 10X 
FCRMAT(1CX 
F C R M A T ( 10 X 
FCRMATC 10X 
F OR V A T ( 1 0 X 
F CRMAT ( OX , 

FOR MAT ( F 1 C 
VR I TF <C . 3 1 7 ) 

WRITE (f. .301 ) 

WRITE <f » 3 C 2 ) 

RFAD(O.lOt) K r 
0 C 10 I = 1 . K 5 

RE AD (0. 1C 1) CL ( I ) .7 I ME ( I ) 


.4 ) 
.( ) 
f> ) 


v; 

K> *.9X,* TIME IN SECONDS •) 


of OCR AH. DATA*,//) 

Hf A T FXCFIANGER - HEAT 


K J 


• r A T f IN 
ex.* MAIN 
OX.* DUAL 

• 1 

maximum t i hr increment, sec 

WINIM'M TIME INCREMENT, SEC 

maximum transfcr fluid flow. 

CX CYCLF D/ T A * .// 3 
f L AO Nf :> TIME. HP 
LOST COLLECTOR HO AT, 
COLLECTOR OUTPUT, KJ 
STORAGE OtTPUT.KJ 
STORAGE INPUT , KJ 
AUX IL I AR V I NPUT • K J 
Excrss heating, kj 

FOILFR RCSUIRFMINT. 

NUMPF R OF TIME INCREMENTS 
CrrKlT RtOUlREvFNT. KJ 
ALLOWABLE POWER FLUID OUTLET 
4 ) 


. F IC, 

, F 1 0 . 

.F 10. 

. F 1 0 . 6 ) 

» F 1 0 • Z ) 

. F 1 0 . 2 ) 

.F 10.2 ) 

.F 10.2 ) 

♦ F 1 0 . E ) 

. F 1 0 . ? ) 

. F 1 0 . 2 ) 

.F 10.2 ) 

» F 1 0 . 2 ) 

. F 1 0 . 4 ) 

.F 10.0 ) 

* 7X » 1 3 ,///) 
. F 1 0 • 2 ) 

TIME • ) 


OF FUSION THERMAL STORA 


¥ V> 


F G/SE C 

. F 1 0 • 1 ) 
.E 10.7 ) 
.LIE. 7 ) 

. E 1 5 . 7 ) 
.E 15.7 ) 
.E 15.7 ) 

. E 1 0 . 7 ) 

. E 1 5 . 7 ) 

,110) 

»E 1 5.7 ) 
TEMP. .K 


• . F10.2 ) 
* .F 10.2 ) 
* . F 10. 4 ) 


• . F10.2 ) 


1C W R 1 T C ( fc . 2 0 1 ) 
TO = T l»‘t (KF ) 

R F A D ( r . 1 0 2 ) 

R. F AD < * . 1 C 3 ) 

R t A 0 ( t , 1 C f, ) 
f f A D ( E . 1 C 6 ) 

W F. 1 T E ( C .314 ) 


CL(I ) . T 1 l, f ( i ) 

CL.CS.CT .CV.TSAT 
TLI N « T R'I f, . 7 PDF S . T AL WPL 
DTm A > ,DT«. I*;, WMAX 
T 0 » T 0 
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P E A D ( E . 1 0 4 ) \» P . 7 1 t A S 7 01- , C 7 K . XN T .ND.PCT 

VPITF ( *> • 2 C 2 ) CL 

WM IF (< ,20 7 ) 7 L I N 

V*PITE(6,?C3) C S 

WPITE{f*2C4) OT 

VPITL (6.2CS) CV 

KPJTE (C , 2 C 6 ) TSAT 

HlTE(6i2C8) 7 P I K 

wriTC (f ,?C9 J TPDF r ' 

V P I T F ( t « 3 3 9 ) TAL^PL 
VP1TF (f .2 10) K P 
kPJTC (t,2 11 ) T 1 
WPITC ( 0 » 3 3 C ) T5 
VPITL~ (0,331) TC 
WPITE(tt213) CTK 
VPITF (f,,2 14 } XNF 

VP ITC (0.212 ) AS TOP 
VRJ7E <C .2 17 ) DTf'IN 
VF1TE(C,216> DTPAX 
VP1TE(5,219) WWAX 
V P 1 7E ( C ,2 19 ) PC T 
VPITCLf .? 15) ND 
F 7 = 0 . 

Orc-WE F=WP* (C TSA7-7P IN )+07+CV-M TPtES-TSA T ) ) 

JC = 1 
¥ 1 = 1 

?e> C NDT=C 
NC= 1 

7L=0.06«ASTGF*CTK/(XNC$CL) 

?PS = 0.C&*AS7 0P*C7K/(XKr*CS ) 

2 PV= 3. Oft* ASTOP«CTK/( XNFSCV ) 

IFCPE.&E.ZPVl CO 70 11 

7 F =ZPV 
CO 70 12 

11 7 P=Z PS 

12 CL CSS =C . 

0 CHG = 0 • 

Cr>7RA=G. 

C S O L A P = 0 • 

C t ORN=C . 

CC ISCH=0 . 

V L =0 • 

L =C 

P 07=0 
K 7= 1 
¥ 2 = 1 

07 ipr=zp/v,p 
C7=D7 I *■' E 

F 7 = E 7407 1VL/3C0C . 

5 C L =L 4 1 

NDT = ND74 1 

CALL FLANAL ( V»L tWF* ,D7 IM( .<1*K7,J6.T1.TD.K2) 

L = C 

El I r (7P0UT • t F.7SA7 ) CO 70 C 2 

rOD!S = »‘P*DTIPr*(CV«(7POUT-TSAT)4CiT4CSA{TSAT-7PIN)) 
C C 7 0 0 0 

E< 2 1 F ( 7 P0U7 .L7.7SA7) CO 7C *3 

CCL IS = WP4DTI Pf «( xnu7^0T+CSn- { 7S AT-7P1 N ) ) 

CO 70 00 

t 3 POP IS = WP->r-71 PF$CS*( 7POUT-7P1N) 

CO cr ISC P=QD ISC H4DU I G 

1 F ( 7PCUT . CT . TPFf S ) OF TO (> 1 
CrUPN = C'L'UPN4CPOV[ ON)! I'-E-DOr is 
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I F { T P0U7 .CT .7ALWFL) GO TO 62 

crfFC7 = QDrFC7-UTM. w‘’L-TPOUT)«!«F $C V«DT I WE 

CO TO 12 

t 1 CTX T PA=GE X T F A+DGO I 5- CPCWf R?DTI ME 
t 2 IF (t< 7 • E G • 1 ) GO TO 7? 

rCCHG = M-»CL*(TLIG-TLOUT)*rTI*'a 
CCPG = GCHG-»DOCHG 
C. SOLAP^GSOLAP + GC^nTI M r 
GLCSO = OLCi S + GC«rTJ ML -DC CMC. 

7 C IF(CT.GE.Tl) GO 7 0 HC 
D 7 1 ME =7 P/fcP 

IF (07 IME.L7 .TTMIN) 0 7 I "C =D T M 1 (si 
IF (DTIME.LE.OTMAY) GC TO 71 
r TlKfc=DT^AX 
71 C 7 = C T + D T IMF 

IF(CT.LE.71 ) GO 7n 7 f 
07 I^ErDTl ME- (CT -T 1 ) 

C 7 =T 1 
* 2-2 

70 ET = ET+DTINE/3GGC . 

K 7= 1 

GO 70 00 

f'u IF (CT.GE.7D) GC Tn 2 0 
IFOC7.L3. 1 ) GC TO A~>c 

00 2 1-2. FO 

I r (CT .GC . 7 I ML ( I ) ) GO TC 2 

CC=QL( 1-1 ) + ( OL ( I )-GL ( 1 -1 ) ) =>(CT-T1ME( I -1 I )/( 7 IMt( I )-T 1ME( I -1 ) ) 

GC TO fc'5 
2 CONTINUE 
Of 0 71 MP=ZP/FP 

1 r(»C7.rQ. 2) W ( =GC/(CL* ( TL TO-TAPPR ) ) 

IF (K7.FQ.?) V>L=CC/(CL=C( TLIN-TLOUT ) ) 

(C 7 = 3 

IF (WL.LE.WMAX) GO TO (7 
VL=W«AX 

°7 DTIML = 7L/V.L 

1 MDT1ML .LT .DT IMP) GC TO 90 
07 IME =D7 1 MP 
GC TO Fe 
9C CTI^c^TIML 
CO 7 0 (8 
ACC C7=T1+DTMAX 

E7 = E7-»DTMAX/3GCC . 

C71Mf=DTMAX 

V 7 = 2 

cc to ro 

f-2 IF (DT IME .LT . 07^ IN) D 7 I ‘■’E =0 7 M I N 
IF (DT IME .LE .L 7MAX ) GC 70 91 
071ME=CTMAX 

91 IF(CT.GE.TC) GO 70 92 
CT=CT4pT IME 

IF(CT.LE.TO) GO 70 A1C 
[ T I ME =[)TIFE-(CT-T r >) 

C 7 = 7 F 

A1C FT=L 74DTIME/3CC0. 

CO TO FO 

92 IF(CT.GE.T6) GO TO 93 
CT=CT4DTI ME 

1F(CT.LE.T6) GO TO A 2 3 
15 TIME =DTI ME - (CT -10 ) 

C 7 = T 6 

A2C r 7=ET-*DTIM1 /36CC. 

GC 70 LO 
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9 3 


43 C 
2 C 


99 


360 


lF(CT.GE.TD) GO TO 20 
C T=CT+DT I ME 

IF(CT.LE.TD) GO TO 430 
DT1ME=DT Jme- (CT-TD ) 
CT=TD 


L T=ET-»DTIML/3COO. 

CO TO SO 
HITE <f .320) 

WRITE (0.32 1 ) CT 
WRITC(6.323) GLOSS 
WRITE (6*324 ) CSCiLAP 
WRITE (6.325) CD1SCH 
WRITE (6.326 ) CCHG 
WRITH6.327) GPURN 
VR ITC (6.336) CDEFCT 
VRjTf (6,326) CEXTRA 
WRITE (6.329) GPCWFR 
WRITE (6.335 ) NDT 
NC = NC-» 1 


IF(NC.GT.ND) GO TO 99 


l=Q 

r co i Fr = *ns( 1 cc . - ioc.«gcmg/qdisch) 
IF(PC01FF .LC.RCT) C,0 TL 9 9 
GO TO 12 

Rl AD ( 5.340 ) AST OR 
JFiASTOR.FQ.l > GO TO ’SO 
J6=2 


CO TO 360 


STOP 
f M> 


EOF ROOT INF ELANAL ( WL.WP.DTIME.K1.K7.J6.T1.TD.K2) 

COMMON TPCUT .TLOOT , T P I N . TL I N , T S AT • OT . CL • CS • CV . CT . XNE . A STD R . XOUT * T A 
#PPR .TF.T6.NC 

CCMMOF/PROPTK/FLM.CPST , STM M, ST ME 

CCMMOF./CALCAO/UAEC ,UAHC »UA HE" .CH.CC .A 1 ,A2 .A3.131 .B2.63.C2 
COMM ON/DI ATRl/Al(*>0)» , 3I(r. C).CI(50)»FI(50)»X1(5C),CY 
COMMON/TMP/ALFA (50 > .Eit TA( SO l.GAMW A (50 ) 

DIMENSION ENER (26) . RSOL L ( 26 ) . T P C 26 ) . XP ( 26 ) , QC ( 26 > 

DIMENSION TL (2C ) ,T{26 ) ,2( 26) 

DIMENSION TEL (2&) .OLOGEC(25) 

10 1 FCf MAT ( 3F 1C . 4 ) 

10? F CEMAT(5F 10. < ) 

103 FORMAT (F 10.2, II »FIO.£ ) 

106 FCRMAT(12X.* SUOROUTINL fLANAL DATA*./) 

IBP FORMAT ( 1 HI .///, • STARTING ELEMENT T t MP • . K *.F1G.2) 

109 F OR M A T ( * STARTING SOLID FRACTION • * F 1 C . 2 ) 

195 FORMAT ( * STARTING TOTAL ENERGY, KJ *.E15.7.///) 

2 4 C F ORMAT(10X,F10.1»10X.F l0.1.1OX,F10.5.10X.F10.2»I0X»F 10.2) 

256 F CRMAT( T1 fc . *E T *. T32 .* CHARGE T 1 ME • . T 52 . • T -66 FLOW », T7 1 . ‘TOLUENE OUT 
(.* ,TV2.*T-66 OUT • ,/) 

4GC FORMAT!* (1) T(I) Z(I) TL ( I ) TP(I) X 

*P(1) *.//) 

46 E F CR M A T ( 1 X . I 3 . 1 X » Fo . 2 . 3 X , F9 • 2 • 3 X » F 9 • 2 • 3 X » F 9 • 2 , 3 X » F 9 . 2 ) 

46 9 FORMAT ( IX * I3.2 5X.F9.2.3X.F9.2.3X.F9.2) 

4 90 F CK MATdX.13.lX.ro. 2, 3X.FO. 2, 15X.F9. 2, 3X.F9. 2) 

49 1 F CR M A T ( IX .I3.37X.F9.2.3X.F9.2) 

492 F OR MAT(1X.* TOTAL SYSTEM ENERGY, KJ *.E15.7,///) 

5 00 FORMAT!* POWER FLUID QUALITY CHECK 500 - Nri= , .I2*» NEF=*.l2) 

501 FORMAT! * POWER FLUID DUALITY CHECK - 501*) 

5 1 C FORMA T ( / / / * • PCWtP FLUID VAPORIZED IN ONE ELEMENT - NO. 510*) 
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fie fcfmat(?il.2E10.a ) 

519 FCFMAT(///»* HEATING FLUID GETTING HOTTER AS IT GOES THROUGH UNIT* 

* ) 

52 0 FORMAT!///.* POWER p LU TO GETTING COOLER AS IT GOES THROUGH UNIT*) 

CS1 FCRMAT(* INTERMEDIATE M4P AT TIME TD= * • F 4 . 1 * • HRS . • ,// ) 

SEE FCRWATt* INTERMEDI ATF PAP AT TIME T 5 = * . F4 . 1 . * H RS . • , // ) 

556 FTRMAT(* INlTrRMf DTATt CAP AT TIME T6 = • , F A . 1 • • HRS . * . // ) 

7C2 FCRMAT(/,* CYCLE FNOING TOTAL ENERGY. KJ *.ri5.7> 

CY = CT 
C 2 = C T 

CC TU(1.?).K1 

1 IF (J6.E0.2) GO TC) 38 C 
N F = XNF 

tl(ne+ i ) = tl in 

Ff=l 

E K I )=-l . 

A I (2 ) = 0. 

C 1(2 )=-l . 

CHGTM=C . 

C T 1 = 0 • 

HF75=T5/3600 . 

H R T 6 = T C / 3 6 0 G • 

FEAD(E.lOl) AFP 1 AFL 1 APL 
RE A D( F . I 02 ) CLFT .XLFFF.XN7.HTCS.HTCV 
FEADEE . 1 C 3 ) PHOC.FG.TIMT 
WRITE (6.6CC ) ArF.AEL * AFL 

GOO FCFMAT(* Arr= ».F10.4,* ael= *.F10.4,» APL= *.F10.4J 
WR1TE(6.6C1) DtFF.XLEFF.XNT 

601 T CRM A T ( • DfFFr *»F1C.6.» XLEFF= *,F10.3.* XNT= *.FI0.3) 
WRITE (6.6C2 ) FHQM . HT C S . HT C V 

602 FORMAT!* RHOV- '.F1C.3.* HTCS= » ,F 1 0 . 6 , • HTCV= • .F 10 .6 ) 

36 C A F PE = AE P/XNE 

AELE=AEL/XNE 
A PLE =A PL /XNE 
AFF E=AEPE -» APLF 
X LE = XLE FF /XNE 
r LM=AS TOR/XNE 
CRXTL = 6 .28 34XLF *XNT 
CPF T=I. 

STFM=C. 

S THE = 0 • 

r=o. 

CALL TKPR CP( TIN I T .EUT L.FPAC. 1 . 1 .2. Jfc ) 

WR1TE(6,188) TIF.IT 
V FITE (6.189) FRAC 
TO 5 1=1. NE 

F NER ( I ) = E NEL 
T ( 1 )=T1N1 T 
r =C+ENER( I ) 

7 ( i ) = F R A C 
5 CONTINUE 
K 1=2 

WRITE (6. 195) E 
WRITE (6.2* 8 ) 

NE 1 =NE A 1 
NLF-NC* 1 

2 KTCL=HL(WL) 

UAHE =HTCL*AELF 
CH=WL#CL 

IF (K6.FQ.2) GO TO 15 

CCS=WP*CS 

CCV=WP#CV 

IF (X 8.00. 1 ) X 6=2 
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IE IF (kL.tO.C. ) GO TO 72 

CUA HC 5= A PL E/ ( 1 ./HTCS+ 1 ./HTCL ) 

CUAHCV=APLE/( 1 ./HTCV+1 ./HTCL ) 

72 EC 3 1=1. NC 

IF (7 ( I ) .EC .0 . ) GO TO 4 

VST A =Z< 1 ) *CLH/(RHOM >XNT ) 

DSTK=SGRT < 1 . 273#VSTK/XLE+DFFF 9 $2 ) 

PEC’LL (I )=ALOG(DSTF/DLFF)/(THCS(T( 1 ))*CR<TL) 

GC TO 3 
4 RSCLL (I ) =C • 

3 C CN T 1 MUE 

GO TO ( 10 . 1C . 1 1 ) ,K7 

F OVER fluid only 

1C TP ( J ) = T F IN 
M4 = 0 

DC 2 0 J= 1 .NE 

>r ( j+ j ) = c . 

IFCTP(J) .GE.TSAT ) GO TO 21 
FLE = RSOLL(J)-»l./(HTCS«APFL ) 
>FTU=1./(RES*WP*CS) 

IF (XNTU.GT . 1 1.0 ) F*FF = . 09900 
IF (XNTL'.GE . 1 1 .5 ) GO TO 22 
E FF=1 .-E XP( -XNTU ) 

22 TP(J+I)=TP(J)+EFF*(T{J)-TP(J)) 

IF(TP{J+1 I.LT.TSAT ) GO TO 1 9 
ETSAT=(TSAT-TP(J))/(TCJ)-TP(J)) 

YNTU=XNTU+ A LOG { 1 TSAI ) 

F TP = YNTU/XNTU 

RtSTP=(RSOLL(J)+l./(HT(XP<J))*APFE))/FTP 

n r i = j 

F XI = 1 .-( YKTU/XNTU) 

GTP=(T(J)-TSAT)/RFSTP 
XP( J+l )=QTP/(WP*QT ) 

I F ( XP { J + 1 ) .LT . 1 . ) GO TO IP 

GTP=Wr»3T 

2 NT U = YNTU £ ( 1 .-1 ,/XP( J + l) ) 

FSUP=7NTU/XNTU 

PESSUr=(RSOLL(J)+l./(HTCV«APFt) l/FSUP 
ZNTU=1 ./ (RESSUP*WP>CV ) 

IF (ZNTU.GF. 1 1.5 ) FFF=. 99999 
IF (ZNTU.GE. 1 l.E>) GO TO 2? 

EFF= 1 . -E X P ( -2NTU) 

23 TI < J+l ) =TSAT+rFF*f T ( J) -TSAT ) 

X F (J+l )=1 • 

CP=QTP + KP< { CV*< TP( J+ 1 ) -TSAT ) +CS*< TSA T-TP( J) ) ) 
F XF = 1 • -F S UP 
NFF = J 
GO TO 35 

18 GP = GT P+CC5=S> ( TSAT-TP( J ) ) 

TP( J+l )=TEAT 

GO TO 35 

19 OF = CC£<"(TF( J+l )*TP( J) ) 

GC TO 35 

21 IF(TP{J) .GT.TSAT) GO TC 25 

Rr.S = RSOLL (J) + l./(HT(XP(J))4AFFr) 

CP=(T( Jl-TSAT J/RES 

XT( J+ 1 )=GF/( WP*CT ) +XP ( J) 

IF < XP( J+ 1 l.LE.l . ) TP { J + 1 )=T SAT 
IF ( XP( J+ I ) ,LE . 1 . ) GO TO 3E 
CF ART = < 1 . -XP ( J) ) AWPOCT 
F >F =GF APT /GP 
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P F 5SUP=< PSOLL ( J Ml . / ( MlCV^APF £ ) )/( 1 . -FXF ) 

N E F = J 

X NTUSF - 1 . /( RE SS UP* WP* C V ) 

IF (XNTU5P.GC .100. ) !' F F = 1 . 

1 F ( XNTUSP.GE .ICO. ) GC 10 20 
t FF = 1 .-CXP( -XNTUGP ) 

26 TF'(J41) = TSAT4CFF#( T ( JJ-T'AT > 

xr ( j+ i )-i . 

C P = G P A P T 4 W P # C V ❖ ( Tr> ( J 4 1 MTS AT ) 

CO TO 35 

25 pfc S=PSOL L (J)4l./(HTCV*APFE) 

>NTU = 1 ./( PESACCV) 

IF (XNTU.GF . IOC. ) F F F = 1 « 

IF (XNTU.Gf . ICC. ) GO TO 2P 
T F F= 1 .-Exr(-XNTU) 

2 5 TMJ+l> = Tr(J)4rFF*(T(J)-TP(J}) 

Cf = CCV*( TP( J+ 1 ) -TP ( J ) ) 

yp( j+i ) = i . 

3 r E NCR (J) = F’NER(J) -OP vD T I M F 
CALL TKPPOP(T{ J) ,TNEP( J) .7( J), 1 .2.2. J6) 

2 C CONTINUE 
PXF = 1 ,-FXF 
CC TO 230 

11 IF (NCI .NE . (NE -* 1 > .OR.NL F .NT . (NE + 1 ) ) GO TO 12 
VP 1TE (6.500) Fri.NEF 
STOP 

PC G I ON I 

NTF IS TOP ELEMENT 

12 CONTINUE 

BOO I F (Nrr .NL .NE ) GO TO 56 

UAEC=PXF/(RSOLL(NT M 1 ./(HTCV*AEPC ) ) 
UAHC=PXF*CUAHCV 
UAHC=PXF#UAHE 
CC=CCV 

CALL AECALCtND 

7LF=A1*T(NE )+A2«TSAT+A3=>lLlN 
TP(nT 4 1 ) = p 1 *T (NE MB2* T 5AT + B3*TL IN 
CCPXFrVP^CV*(TP(Nir4 1 MTS AT ) 

CC TO 55 
C 

C NCF IS NCXT-TO-TQP ELEMENT 

56 IF (NEF.NC . ( NC- 1 ) ) CO TC Ctt 
U A HE = P' X F * U A H F 
UAHC = PXF>>CUAHCV 

UAEC = RXF/(PSCLL ( NE-1 )+ 1 ./ ( HTCV« AE PL ) ) 

CC=CCV 

CALL APCALC ( NE - 1 ) 

Cl ( 1 ) = A 3 

F I ( 1 ) =-A 1 *T (NC- 1 ) -A?*TSA T 
PI (2 ) = B 3 

F 1 ( 2 ) = -fc 1 *T (NT - 1 ) -l' 2* T' AT 
U A HE =HTCL*ACLf 
UAHC=-CUfiHCV 

UAFC=l./( RSOl L ( f.F ) + 1 ./(HTCV*Afl E ) ) 

CALL AFCALC(NC) 

A I (3 >=-l . 
f 1 (3 )=A2 
C 1 ( 3 ) =0. 

f 1 (3 )=-A 1 < T ( NC )-A3* 1 L IN 
A I ( A ) = 0 2 
P I ( A ) =- 1 . 
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on 


F 1(4 )=-Bl*T(NT )-F7*TL IN 
CALL TP ID 1A ( 4 ) 

TLF = X 1(1) 

TL(NE)=X1 (2) 

TMNE ) = X I ( 3 ) 

7 F ( NF 4 1 ) = XI ( 4 ) 
CCFXFrKP*CV*(TP(NF)-TSAT > 

CCINT )=WP*CV*(TP(NE4 1 ) -7P(NE ) ) 
CO 70 55 


RTGICN I IS 3 OR MORE CLEMENTS 
56 UAHC=PXF*CUAhCV 

L) AHE =FXF*HTCL*AE LF 

UAT C-PXF/(RSOLL (NFF ) 4 j ./{HTCV*ACPE ) ) 

C C = CCV 

CALL AFCALC(NEF) 

C I ( 1 ) = A3 

F 1 ( 1)=-A 1*T (KEF ) - 42->7SAT 

r 1 { 2 ) = [ 3 

FI(2)=-D1*T(NCF}-P2*7SA7 

UAHC=CUAHCV 

UAHE=H7CL*AFLr 

L 1 = NF F4 1 

L2=NE - 1 

J = 1 

DO 50 I=L1.L2 

UAEC = 1./(RS0LL( I )4 l./(HTCV«AFPE )) 

CALL AFCALC ( I ) 

J = J4 J 

A ](2*J-1 ) = -l . 

E I(2*J-1 )=A2 
CI(2*J-1 ) = A 3 
r I (2* J- 1 }=-A 1*T( 1 ) 

A 1(2*J)=B2 
PI (2*J)=B3 
CI(2*J)=-1. 

50 F I (2<«J) = -E 1*7(1) 

UAEC=1 ./( RSOLL (NT ) 41 ,/(H7CV*AEF f ) ) 

CALL ALCALC( NT ) 

J = J4 l 

A 1 (2* J-l )=-l . 

F I (2* J- 1 )=A2 
C 1 ( 2* J- 1 ) =0 • 

F I (2* J-l )=-A 1*T(NF )-A 3*7LIN 
A I( 2* J)=P2 
P 1 ( 2 * J ) = - 1 . 

r 1(2*J) = -L l*7(Nt )-B3*7LIN 
N = ? * ( NE -NF F 4 1 ) 

CALL TRIDIA(N) 

7LF = X I ( 1 ) 

7L(NLF4J )=XI (2) 

J = 1 

7 F (KEF ) = 7 SAT 
1. l=Nr F4 1 
L2=NE-1 
DO 5 1 1 = L 1 t L 2 

J = J + 1 

TF( 1 )=XI (2*J-1 ) 

CC(I-l)-HP*CV*(TP(l)-Tr(I-J )) 
TL(I+1)=XI (2*J) 

51 CTNTINUE 

7 f ( N'E (-X 1 (N - 1 ) 

7F(NL41)-XI (N) 
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GC (NE.-l )=WP=>CV*(TP(NE ) -TD(NE“1 ) ) 
CC(NE ) = WF *CV»(TP(NF+1 )-TP(N'E ) ) 
OCFXF=QC(NEF ) 


REGION II 


SECTION FXF 

ti IF(AXF.EG.l) GO TO r >9 
FXr = 1 • -R X F 

CUAHC7 = APLF/ ( 1 ./hT (XP(NL'F ) ) 4 l./HTCL) 

UAHCF = F XF SCUAHCT 
UAHEF=FXF$HTCL*AELC 

OAF cr=FXF/(RSCLL(rrF ) + l»/(HT (XP(NCF) )*AEPE ) ) 

CALL VAP (UAECr .UAMC * 7 . U AHE F , C H, T S A T . T ( NETF ) * T L F , OHr . GE C . GHC . T L ( NE F ) ) 

OLCGLC(NEF)=crC 

f. CtNEF ) = 3HC+CrC + GCRXF 

xr( NEF )= I . - ( CHC4 0FC) / ( WP<-CT ) 

If(XP(NCF ) .GE .0. ) GO TO SR 


NFF IS NCI GOING DOWN 

F =C 

>I1=G. 

X 1 «- — 1 • 

>Pl=XP(NrF ) 

XF 2 = 1 . 

It t<=X+l 

IF (X. EC. 3) GC TO 7<* 

X1? = XI 1-XP1* (XI2-X 1 1 )/(XP?-XPl ) 

F X 1 = ( l.-XI3)<TXF 
UAHC 1 =RX I *CUAHC T 
UAHLI=RXI»HTCL*A[ LE 

UAFC I=RX1 /(RSOLL CNLF )4 l./(HT (XP (NEF) )*AEPE ) ) 

C ALL VAP( UAE Cl .UAHCl ,UAHF 1 ,CH. T SAT *T ( NEF ) • T LF • CHE .GE C.GHC.TL I ) 

CLDGEC (NEF) =CCC 

CCRX I =0E C4GHC4QCRXF 

XPX I 3= 1.- (0LC4GHC )/( W'P*OT ) 

IF (XP> 13) 76.7P.77 

7 e XF1=XFXI3 
X 1 1 =XI3 
CC TO 70 
77 > P£ = X FX 13 
X I2 = X13 
C C TO 7S 
76 FXI = X 1?$FXF 
N F I = N F F 
1 =NLF 
GC TO 79 

RTGIGN II MIDDLE ELEMENTS 
IQ IF(NEF.EQ.l) GO TO 270 
1 =NEF 

L‘AhC2 = CUAHCT 
UAHE? = HTCL=>AELE 
60 1=1-1 

UAFC? = 1 . / ( R S DLL ( I ) 4 I . / < M T ( XP ( I ) ) * A E PE ) ) 

CALL VAP( UAEC2 .UAMC2 ,U AHT? .CH.TSAT ,T ( I ) .TL ( 14 1 ) » QHE . GE C * QHC * TL ( I ) ) 
Cl OGEC( 1 )=QEC 
GC ( 1 )=GHC4QEC 

XP(1 )=XP(1+1 )-(GMC4'JC C)/(WP*GT ) 

T T ( 1 ) =T S A T 

1F( XP (I ) .GT .XP( 14 1 ) ) GC TO ei 
1 F ( I . FQ. 1 . ANr .XF ( 1 ). C,F .0. ) GO TO 0 1 
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if (xr< i ) . gt . 0. ) on tc go 

REGlCN lit IN NCI. FIND RXI 

X =C 

X I I =0 . 

X 12 = 1 . 

X P 1 = XF ( I ) 

X P2 = XT ( 1 4 l ) 

7 C K=K+1 

1F(K.EG.3) GC TD *0 
X13=XI1-XP1*(XI2-XI1 )/<XP?-XPl ) 

PXI=1 .-XI 7 
UAHCI =RXl SCUAMCT 
UAHCI=RX l*HTCL«ALl.L 

U A I C I =RX 1 / ( RSCLL ( I ) + 1 »/(HT (XP( 1 ) ) # AE PE ) ) 

CALL VAP( UAEC I .UAHCI , U AME I , CM , TSA T • T ( I ).TL( Ml ). GhE . Or C .GHC . TLI ) 
CLC GE C<NE I ) = GEC 
CCRX I=GFC4GHC 

XPXI 3=XP( 1 + 1 )-( GFC4GHC )/{ Vt>#QT ) 

IF (XPXI3) 71 .or ,73 
7 1 X P 1 = X PX 1 3 
X I1=X 13 
CL TO 70 
73 XP2=XPXI3 
X I? = X I 3 
60 TO 70 
ec F X I = X 13 
N E I = I 
CC TO 7y 


SF T CONDITION Ml 

FI L3=I+1 
V 1 = 1 
F X I = 1 . 

CLDGEC (ME I ) =C . 

CCRXI=0. 

T L I = TL ( I-U) 

NE 1 =1 
XP( I ) = 0 • 

60 TO 75 


270 1=1 

TL I=TLF 
M 3= 1 

CCR XI =0CRXF 
FXI = 1 ,-RXF 
NE I=NTF 


ST T CONDITION M3 


REGION III 

REGION III ID CNF ELEMENT 
79 I F ( 1 .NE . 1 ) GO TC 03 
L'AMC=FXI C CUAHC5 
UAHC=FX I <• HT CL^AELF 

UAEC =F XI /(RSCLL < 1 Ml ♦/(MTCS’XAERD) 

CC=CCS 

CALL ADCALC ( 1 ) 

T L ( 1 ) = A 1 * T ( 1 ) 4 a P i-r P I N4 a 3 => T L T 
TTR I AL=B1 *T < 1 >4P?*TP IN + H 3-7 TL I 
OC ( 1 ) = Wf «CS« ( TTR I AL- TP IN )4QCRX1 
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IF CH.CG.l ) TP(2)=TTPIAL 
T F ( 1 ) = 7 P I N 
C,C TO 15C 
C 

C REGION III IS Tun LLCMINTS 

PC IF (I. Nr. 2) GT TO 100 
UAHC=CUAHCS 
U A HE =HTCL*A£LE 

UAEC=1./(RS01L( 1 )•»!./( HTCS*ACPE )) 

CC=CCS 

CALL A F C A L C ( 1 ) 

C 1 ( 1 )=A 3 

F I( 1 )=-A 1 *T( 1 ) - A 2 * T P I N 
E I < 2 ) =B3 
F I (2 ) =-B 1 *T < 1 
UAHC=FX1*CUAHCS 
UAHC=FXI*HTCL*AELF 

UAL C=FX I /< PS DLL ( ? ) ♦ 1 ./ (H TC S*AC F C > ) 

CALL APCALC( 2 ) 

A I ( 3 ) = -l . 
r> i ( o ) =a2 
C I (3 )=C . 

F I <3 ) =-A 1 *T ( 2 )-A 3*TL I 
AJ(A) =P2 
r I (A )=-l . 

F I ( A ) =-B 1 *T ( 2 )-EO*TLI 
CALL TRIDIA(A) 

TL ( 1 >=> 1 ( 1) 

TL(2)=XI (2) 

Tr <2 > I ( 3) 

TTR I A L = X I (A ) 

IF (Ml .LG. 1 ) TP( 3 )=TTP I AL 
OC { 1 )=WP*CS* (TP( 2 ) -TP I N) 

CC( 2 )=VP*CS*< TTR IAL-TP (2 > HCCPX I 

CO TO ISC 

REGION III IS T H P T E OP MOPE ELEMENTS 
IOC UAHC=CUAHCS 

UAHE=H7CL AAELL 

UAEC = 1 ./< PS DLL ( 1 )* 1 ./ (HTCS^ATPC ) ) 

CC=CCS 

CALL A PC ALC ( 1 ) 

Cl ( 1 )=A3 

F I ( 1 ) --A 1 *T ( 1 )-A?*TPlN 
t 1(2)=F 3 

F1(2)=-P1*T( 1 )-F'2^TPIM 
L2 = NF 1-1 
OC 110 1 = 2 . L 2 

UAf C = 1 ./(PSOLL( I )■*■!. /(HTCS*AFPt )) 

CALL AECALC ( 1 ) 

AI(2*I-1)=-1 . 
f I ( 2 * 1 - 1 ) = A 2 
CK2AI-1 ) = A 3 
F I (2*1-1 ) =- A 1*T< 1 ) 

A ] (2*1 ) = P2 
f I (2*1 )=P3 
CI(2*I) = - 1 . 
ii o Fi(?*n = -n*T(i) 

U A H C = F X 1 *CUAHCS 
UAHL=FXl*HTCL*AFL r 

UA f C = FX1/(P SOL L ( NE I ) ♦ 1 ./ (HTCS*ALPT )) 

CALL AECALC ( NE I ) 

A 1 (?<NE I - 1 > = -l . 


D-2 3 



c 

c 

c 

c 


c 

c 


c 

c 


c 

c 


l 1 ( 2$FvE I - 1 ) = A f' 

C1(2$NE 1- 1 >=C. 

fi( 2 *kfi-i)=-ai«t(nf;i)-^?*tli 

A 1 (2*NE 1 ) = B2 
r i (2*Nr l ) = - i . 

r I (2*NE I ) = -B 1*T ( NT I ) -0?$TL I 
CALL 7R1DIA(NE 1*2) 

TL( 1>=X I ( n 

TL ( 2 )=XI ( 2 ) 

L2 = NF I- 1 
TC 12C 1=2. L2 
TK I ) = XI ( 2*1 -1 ) 

CC(I-1) = *F*CS$<7P<I)-TMI-1)> 

120 TL(1+1)=X1(2«1) 

TF(NF 1 )=XI(2*N£I-1 ) 

T7F IAL = X1 (2*NL 1 ) 

IF(Ml.LQ.l) 7P(NF 1+ 1 )=TTR1AL 

CC (NF 1 -1 ) = t»P$CF>« ( T P( NF I ) -TP{ NL I -1 ) ) 

r,C(NFl)=WF*CS->(TTPML-TP(NEI))4CCRXl 


IOC CCN7INtr 

REGION III - CODING UP 

P 1 = C 

NF ] 7MP = Kjf I 
NT f TMP=NE F 

1 F (NL 1 .E C . 1 ) GO 71 1 72 
l l = Nfc 1-1 
PC 16 0 1 = 1 » L 2 

tCF = (7L< !-» 1 >-7L (J ) )->V*L<CL-(7P(I + l)-7P(l) )*WP$CS 
XP( 11=0. 

FNFR( I )=ENE R( I ) 4 D 0 F -> D 7 I ME 

CALL 7KPRGP ( 7 ( I ) .'"NCR ( I ) ,7 ( I ) .2 .2 .2. Jfc ) 

IF (7PCI41 ) .G7.7EAT ) CO 70 166 

16 0 CCN7INUE 
CF 70 172 

7 S A 7 OCCURS f F F 0 R E NF I 
If b X P ( J + 1 ) = CS6 (7P( 1+ 1 1-7S/7 )/Q7 
I F ( XP ( I 4 1 ) . G 7 . 1 . ) GO 7 0 16 7 

7P( I + 1 1 = 7 S A T 
NL I7MP=I 
CC 70 1£0 

VAPOR 12A 7 ION OCCURS IN ONl FlFMFN 7 
If 7 NT F 7 M 1 
NFI7MP=I 

PXF = CV6 ( 7P ( I 4 l )- 7 SA 7 )/ (C5* ( 7P( 14 1) -7P( I )> ) 

TM I 4 1 ) = ( CS# ( 7P( 14 l ) - 7 SAT )-Q7 )/CV47SA7 
X P ( 1 4 l ) = l . 

GO TO 2 1 C 

CLTPCNT NO I 

172 I-NEIIMP 

XF<N(:iTMP4I) = (0C<Nri7MP)-WP^CS*(TSA7-TP(NE I TMF )) )/( WPCQ7 ) 

IF ( XP (NE I 7MP4 1 ) ,l.T # ) co TO 17A 

IF (XP(NFITMP+1 ).G T . 1. ) GO TO 176 

7 P ( NC ITMP4 j ) = T 5 A T 

0GF = WL*CLMTL<I4l)-TL(l))-GC(I) 

F NF R { I ) = F M. R ( 1 )4D0L 'H T I Mf 

CALL TKPR OP ( T< I ) .T NEP ( I ) ,2 ( I ) .2 .2 .2 . J6 ) 

GO 7 C IRC 
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IN NL I . If AT DIC NOT OCCUR 

17 A XP(I+ 1 )=C. 

TP( I + 1 )=TP ( I I+OC (NE I T M f ) /( rfPsCS ) 

N r 1 TMF=Nf I TMP + 1 

rcc = v.L «CL* < TL ( I -» 1 )-TL( I ) ) -oc ( 1 ) 

FN[P(I)=ENrR(I)+DnF->rTIvr 

C ALL TtkPPTP ( T ( I ) ,r NE R ( I > ,Z ( I ) . 2.2 .2 . Jt ) 

C.C TO 172 


IN MCI. T CANS IT ION IN ONE F L F- PL NT 

176 TP(NL IT VP-H ) = (OC CNF I TPP) -( TsaT-TP(NE I 1MP ) ) * WPeC S - WP* QT )/ (WP*Cv )+T£ 
1 A T 

l-FT T MP= Nt I IMP 

I F CNLFTvp.E 0 .(IF 1 ) GC. TP ?C1 
350 IF (NEFTMP.EQ.Nrr ) CO TP 303 

P XF = WP*CV* ( TF (NE 1 TMP4 1 )-TSAT ) /OC (NE I TMP ) 

D(,r= ( TL ( 1 + 1 ) -TL( I ) )->»LPCL-QC (NF I ) 

E NF R ( 1)=ENER< I ) +DGL*L'T 1ME 

CALL TNPRLP (TCI ) .FNTR ( 1 ) ,7< I ) .2.2 .2. JC ) 

C C T C, 21C 

NfFTMP IS MI 

3C 1 IFCNEFTPF .EO.NLF) GO TO 3C2 
05 1 CA = GC (NLI ) -QCPX I 

yfA=XP(NFFTMP-H)-OCRXI/(Wr«OT) 

I F ( XPA- 1 . ) 310.311,312 

XP A TGU/Lf ONL 

311 R’XF=1.-FXI 
TPCNFFTMP+l ) - TSAT + (JCPX I/CCV 
GO TO 305 

XPA IS LESS THAN ONE 
310 C0= ( XP< NE F TMF4 1 ) - 1 . ) ❖WF SOT 
P XFrCD/QCPX I <•( 1 .-F XI ) 

TP (NF F TM3 -» 1 ) = TSAT4 OD /CCV 
Cf TO ?0t 

XPA IS G R T A T L R THAN ONE 

312 00 = ( XFA-1 . ) ❖ ViFsCT 
PXF=0C/aA6FXl41 ,-FXI 
GSUPHT = QC4 OCKX 1 

TF (NLFTMP41 )=TSAT-»GSUPHT/CCV 
C C TO 3 C 5 

f.r r TMP IS NF I 6 NE F 
30 2 IF(M3.EG .1) GO TO 3 3 
GP=GT#*P 

G'A=OC (NE F TMP ) -OF - OCRXr 
XPE = XF ( NEF THP4 1 )-r,CRXF/( WP*GT ) 

IF(XP{-1.) 3CA»3(i 5,306 

XPP l r S S THAN ONL 

30 A 0 F Li P H T =• ( XP (NPF TPP4 1)-1 • )-XGT«WP 
F X.F =RXF6GSUPHT/GCPXF 
TKNfFTMP-t 1 ) = TSA T + QSUPHT/CCV 
GC TO 3 0 f 

XPP GRCATFR THAN ONL 
3C ( XPA = XP( Nt FT HP 4 1 ) -OCRX I/( wfP*OT ) 

IF(XPA- 1 .) 307.30°, ->06 
C 
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XPA IE ONT 


C 

c 


c 

c 


3Cf r »f = i . -Fy i 

GSUPHT=GP4QCRXF 
Tr<NEF7MP4 1 ) =TSA74QSUPHT/CCV 
CT TO 2 Cf 

XPA IS LESS THAN ONT 
307 Cl =( XTB- 1 . ) *WP*CT 

P XF = RXF + ( GD/OF )<•( 1 . -RXF-FXI ) 

C SUPH T=GD4QCPXF 

TP(NEFTvP4l )=TSAT+03UPHT/CCV 
XPP IS ONF 

30 6 ret = ( TL <NF TTFP+ 1 )- TL ( NEF TMP ) ) ’> V, L« CL - OC ( N L F T«P ) 
r EER ( NEF TMP ) = ENE P( NT FT^P)-«DQF«DTI ME 

CALL TPPPOPt 7 (NFFTMP ) » LNER ( NETFT MP ) * Z ( NEF TMP ) *2 »2 *2 • 36 ) 
XP( I 4 1 ) = 1 • 

CP TO 210 


XPA IS CRFATER THAN ONE 
300 Ot=(XfA-l .)«WP*CT 

PXF= 1 ( 1 • “ QD/G A 1 *F X I 

C6UPH7 = OD«GB-»GCRXr 
Tr(NrrTMP-H)=TSAT-»OSUPHT/CCV 
CC TO 306 


NF F TMP = NE F . NOT NE I 
303 CA=QC(NEF I-OCRXF 

XF A = XF < NE FT MP4 1 ) -OCGXF/F WP *GT ) 

I F ( X P A - 1 . ) 32C.3C r ..?22 


XPA IS LESS THAN 1 
12 0 Q SUPH T = ( XP (Nf_rTMP4 1 ) - 1 . ISWPSQT 
P XF=RXF#G5UPHT/GCPXF 
TP(NEFTMP41 ) =TSAT40SUPHT/CCV 
CO TO 3C6 


XPA IS GREATEF* THAN 1 
??2 cr= c xp a- i . >$cT*wr 
P XF = PXF +F XF $ GD/C A 
0 SUPH T = GD4QCPXF 

TP(NCFTMP-*1 ) =TSAT 4CSUPHT/CCV 
CC TO 30 6. 


RFG ION II. COVING UP 


100 1=141 

1 F ( I .NE .NL4 1 ) GO TO 101 
W P I T T ( t . 6 G 1 ) 
r top 

101 X F C I 4 1 ) = xr (I )4QC(T )/(*’P*GT) 

PGE = (TL< I 4 l )-7L( 1 > )$V*L*CL”GC( I ) 

1 F ( XP< 1 4 1 ) .GT. 1 . ) GO TO 2CC 

E NER ( 1 )=ENCR { I ) + F 0 E * l> 7 1 M E 

C all TKPROP ( T ( I ) ,F NER < 1 ) . ? ( I ) ,2 .2 .2. JC ) 

GO TO 190 

In r L f v t NT NEF 

2 r, 0 MFTMPrI 

TP(I41)=(GC(l)-VP«CT*{l.-XP(I)))/( WP*CV } 4 t s a t 
IF (NEF TMP .C 3 .NT I ) GU TO 301 
362 If (NEFTMP .FO.N'E F ) GO TC 3C3 
y m 1 4 1 ) = i . 
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rf.KR (I) = ENCR{ 1 MOW «[.T I“E 

CALL TKPRCP( 7 < I ) .FNER < I > .2 < 1 ) • 2.2.2, J6 ) 
RXr=l.-WP*QT*(l.-yP(I) )/QC ( I ) 


REG I PN It COMING UP 


2 16 C T N T 1 NUE 
PP1 L1=NEFTMF-*1 
S UM = 0 . 

rr 155 J = NE I » NF F 

IF (( TFL( I )-T6AT ) .1.7 . IE-1 0) GO 70 155 
SUM=5UV+< l.-(T(I)-TSAT)/(TEL(I)-TSAT) ) SOLDO EC ( 1 ) 
165 CCK7 1MJF 

CANT = (OLDGEC ( NET )*( T ( I >-T5aT >/( TEL <1 )-TSAT > )/FXF 
pPXr=-6UM/0ANT 

rxf = rxf-*drxf 

66 8 IF (RXF.LE ,C. ) N E F T M P = N T F 7 M P + 1 
IF (RXF.LE .0. ) RXF = RXF-»1. 

IF (RXF.LE .0. ) GO 70 665 
665 IF (RXF.6T. 1 . ) NT FTMP=ur F 7MP- 1 
IF (RXF .37 .1 • ) RXF=RXF-1. 

1 F ( RXF. G7 . 1 . ) GC 70 665 
IF (RXF.L7..5E -1 ) PXF=.CE-1 
I F (L 1 .GT . NE > GO 70 229 
DC 220 I = L 1 t NE 

rCF^WLSCLSI 7L( I + 1 ) - T L ( 1 ) ) - GC ( 1 ) 

ENER ( I )=ENER ( I ) +OOEXD7 1ME 

CALL 7ICPR0P ( T ( I ) .ENER ( I ) »Z ( I ) .2 .2 .2. J6 ) 

22 C CONTINUE 

C 

229 CONTINUE 

NF I = NF I IMP 
Nf F=NEF7MP 
236 TPOUT = TP(NE+ 1 ) 

DC 1 5 A K- 1 • NE 
16A TFL(K)=T(K) 

7L0UT = T L ( 1 ) 

N 3 — C 

X OUT = XF ( N F 1 ) 

E7=CT/360C. 

IF (CT.FQ.T5) GO TO 2 C f 
IF (CT.EG.T6) GO TO 206 
IF (CT.EO.TD) GO 7 r 7CC 
C7 1-CT 1 ♦ D7 IMF 

1 F(C7 .GE . T5. AND .C7 .LC . 76 ) CHGTF = ( CT-T5 )/3c 00. 

IF (CT1.LT .1762. ) GO TO 2 r 7 

WRITE (A » 2 A 0 ) r T » CNGT M t WL t 7 POUT t 7 LOUT 

C 7 1 =0 . 

257 RETURN 

25 L WRITE (fc.f r 5 ) HRT 5 
GO TO 26 C 

25C WPITE(6,6F6) HRT 6 

CHCTM= ( Tf -7 5 )/36 C'P . 

2 (A F=C. 

WR] TF ( 6 « A 0 0 ) 

FT 2 5 2 1 - 1 t N E 

IF (XP(I ) . L T . 0 • ) XP( I ) = C. 

WRITF ( t »A( 8 ) I .T( I ) . Z ( I) tTL ( 1 ) t TP ( I ) .XP( 1 ) 

2‘ 2 F =E + ENE R ( 1 ) 

J =NE ♦ 1 

WRITE ( f .AFP) J . T L ( J ) .71 ( J ) • XF ( J ) 

W P 1 T t ( 6 . A 9 2 ) T 
WRITE (6 .26 5 ) 
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vmf(6.2«0) CT.CHCTM.WL.TPOUT.TLOUT 
C HCT M =0 • 

C T 1 = 0 . 

T A PPR = T ( 1 ) 

R E TURN 

7 c o r=c. 

C'C 70 1 1-1 » NE 
70 1 F =E + CNE R ( I ) 

W R I 7 E ( 6 . 7 C 2 ) f 
C 7 1 = 0 . 

RETURN 
E ND 


SUBROUTINE ALCALC(l) 

CCMMON/CALCAD/UAEC . JAHC.UA HE .CH.CC.A1,A2.A3.B1,B2.B3.C2 

THIS SUBROUTINE PROVIDES COEFFICIENTS FOR THE SIMULTANEOUS LINEAR 

ECUAT10NS WHICH SOLVE THE THREE-WAY ANALYSIS. 

ALFA 1 = (UAhE n AHC )/CH-< UAHC + UAE t )/CC 

AIF A 2= ( UAKC**2- (UAHC+U AHE ) * (UAE C*UAHC ) )/( CH*CC ) 

TEST=ALFA1**2-A*AIFA2 

rriKCRl = (- ALFA14S0PT ( Tf ST ) )/2. 

R0WER?= (- ALF A l-SQRT ( Tt ST ) )/2. 

Z 1=P0WER1 *CH/UAHC* 1 • 4 U/ Hf /UA HC 
22=PCWER2*CH/UAHC4 1 . 4UAHE/UAHC 
1 F (POWER 1 .LT.-1AE. ) P0WER1=-1AE. 

1F(P0I*LR2.LT • - 1 A ft • ) P0WCR2=-1AE. 
rrwRi=EXP(powrRi ) 

F PWR2=E XP (P0WER2 ) 

FC IFF =E PWP 1 -LPWP2 
X = A l?S ( POH f R 1 ) 

Y-A6S(PDHFR2) 

IFtX.CT.l ,E-2.AF:C' .Y.GT.I.E-?) CO TO A 
FDIFF=PDWER 1 +POW F P 1 * *2 /2 . +P0 WF.R I **3/6 . 

E riFF=EO I FF-P0WLR?-P0wr R2**2/2 . -P0WER2**3/6 . 
a e-zi*epwri-z?«e twr2 

A 1 = 1 .♦(EE'HRl * C P H R 2 * ( Z2 -Z 1 ) -ED IF F )/B 
A 2=ED IF F /ft 

A 3=EPHRl*rPWR2* ( Z 1 “Z2 ) /B 
Pl=l.-(Zl-Z2+EDlFF*Zl*22)/b 
F2= ( 2 1-22 )/E3 
F'3=2 1*22*EDI FF/8 
R f TURN 
F ND 


SUBROUTINE V AP ( U AF C . U A HC . U A HC . C H. T S A T . TE . T L 1 .QHE.OEC .OHC.TL?) 
F CWCR=< UAHE4UAHC )/CH 
IE(POWER.GT. 1A5. ) P0WfR=145. 

F PWR=EXP< -POWER ) 

CDNST 1= ( UAHE *TF. + UAHC*T SAT )/(UAHE+ UAHC ) 

CDNS T 2= CH* ( TLl-CONSTl ) /(UAHE+UAHC ) 

CHE=UAHE*(C0NST2*( 1 .-E PWP )-» CONS T1 - TE ) 

0 HC - UA HC* ( CONST 2* ( 1 . -F PWR ) 4 CONE T 1 -TS AT ) 

CFC=UAE C* ( TE-TSAT ) 

TL2=(TL1-C0NST1 )*FPWR4C0N57 1 
RFTURN 
E NO 
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SLTR0UT1M TKJDM(N) 

CtMMON/DI ATR1/A1 (50) .F. I( EG ) .Cl (5C ) »F I (50 ) .XII50) *CY 
CCNMON/T^p/ALF A ( f 0 ) « E't T A ( fcC ) . G A MM A ( 50 ) 
lMGI(l).fG.O.) an TO 100 
f<fta( i >=r i ( l ) 

C A WM A(1 1 = C I ( 1 )/( ETA{ 1 ) 

NN1=N-1 
TO 10 1=1 .NM1 

IF(BElA( 1 ).EC.O. ) GO TL 101 
GAMMA (I ) = C I < 1 )/EFTa( 1 ) 

1M = 1«1 

10 FE TA ( IP 1 ) = P 1 ( IP1 )-A 1 ( IP1 )*GAMMA ( I ) 

ALT A ( 1) = F 1 ( 1 )/Pt T A ( 1 ) 

DC 30 I = 2.N 

3 0 ALF A ( I)r (F 1 ( 1)-A I ( I ) *ALF A( I -1 ) )/BE TA ( 1 ) 

X 1 <N ) =ALF A ( N ) 

DCt AC 1= 1 .K w 1 
J = N-I 

AC X 1 ( J ) =ALF A ( J )-GA ^MA( J)*X 1 ( J+l ) 

PE TURN 

ICO VPllFlCif CO) 

CGC FOPMATC »b I ( 1 )=0 IN T R 1 D 1 A * ) 

STCP 

1C1 VF'llE (6.CC1 ) 1 

60 1 FORMATE* BCTA(I) IN TR1D1A=0 FCR 1=».12) 

STCP 
E ND 


SUBROUTINE TtCTRCP ( TEM1 , C N TH , FP A C . M ,K 2 .N . J6 ) 

THIS SUBROUT 1 ME ESTABLISHES TABLES OF CLEMENT CNTHALPY VS. THF 
FRACTION OF SOLID THCRHKEfP AND TEMPERATURE 
COMMON/PROPTK/FLN.CPET .STMM.STME 
DIMENSION E M (1S)»FE( 16).T(18).H(18).Z(18) 

200 FCFMAT(12) 

2C1 FORMAT (F 1C. 2 .F12.B.F 10 .6) 

202 F 0PMA7 ( /,* X . • TAELE OF TE MPFRA T URE .E NTHALPV . AND SOLID FRACTION 
fcOF THCRMKTEP • ./) 

203 F0RMAT(10X.« T(1)*.10X,» H(l)*.lCX,» Z ( I ) * ./ ) 

20 A FORMATE 6X ,F 10 .2.5X ,F 12.5 .EX ,F 11 .6 ) 

2CS F OPM A 7 ( 5X « * ENERGY BEYOND TABLE LIMITS FCR E NT H= • , E 1 0 . 4 ) 

IF (N.ECi.l ) GO TC TOO 
GC 70(1. 2), <2 
1 IF (J6.ro. 2) GO TO 10 1 
VF1TE (6.202) 

READ(E .200) K4 
VP1TE (6.203 ) 

10 1 J4=K4-1 

J3= 1 

DC 11 1 = 1 .K 4 

IF<J6.EU.?) GO TO 102 

PE AD (f. 201) T(I),H(i). 2(1) 

1C2 CONTINUE 

IF(K2.E0«2) GO TO 2 
IF(J3.E0«2) GO TO 10 
IF (T( 1 ) .L7.TEMP) GO TO 10 

ENTH=rLMMHC I -1 )4( TSMP-T ( I -1 ) ) #(H ( I )-H(I -1 I )/( 7 ( I )-T ( I - 1 ) ) ) 

F RAC=Z( I )-K T( 1 ) -TFMP ) v (Z ( I - 1 1 -Z ( I ) )/( T( I )-T ( 1-1 ) ) 

J3=? 

10 EKI ) = E L M <■ H ( 1 ) -*STMM=>CPST<-7 ( I ) 

FE ( I )=ELM*H( 1 > + STME*CR5T*T( 1) 

I F ( J ( . E Q . 2 ) CD TO 11 


D-29 



no non non nnn 


M.lTL(ft.rC4) TC I ) ,H( I ) .2 ( I ) 

11 CONTINUE 
RETURN 

2 Cl T0<3,4 )»M 

IF (l NTH .CT. r M(<* ) .r«.CNTH.LT.EM{ 1 ) ) GO TO 25 C 
DC 0 1=1, JA 

IF (CNTH.GT.EMt Ml ) ) GU TO 5 

T F NP = T ( I ) 4 (T( I + 1 ) ” T ( I ) )*(ENTH-EM( 1 ) ) / ( E M ( I 4- 1 ) - F.M ( 1 ) ) 
FRAC=7(I)4(2(I+1)-2(I))»(L : NTH-EM(I))/(EMCI+1)-EM(I)) 
RETURN 
5 CCNTINUE 
GO TO 2 5 C 

A 1F(ENTH.GT.FE(KA ) • 0 R • E N T M • L T . T E ( 1 ) ) GO TO 250 
DC 6 1=1. JA 

ir<E N7H.GT.EE (1+1)) GO TO 6 

TEMP=T(I>4(T(I + 1)-T(I) IMENTH-EEF I ) )/<CE ( 1 4 1 )-EE ( J ) ) 

FRAC = 7( I >4 (z < 14 1 )-Z < I ) )*(EnTH-EE( 1 ) )/<CE( 14 i ) -EE 1 1 ) > 
RETURN 
f CONTINUE 
GO TO 303 
3C0 DC 302 1= 1 ,KA 

IF (T< I l.LT.TFMD) GO TO 302 

FR/C = Z( I J 4 C T ( I)-TF MP )< (Z(I-1)“Z(1))/(T(I)“T(I“1)) 
IF(M.E0.2) GO TO 711 

E NTH=FM( I-l)4(EM(I)-tM(l-l) )«(TEMP-T <1-1))/(T(I)-T<1-1>) 
GT TO 301 

31 I FNTH = EE ( 1 - 1 ) 4 ( C l ( I ) -C F ( I - 1 ) ) * ( TEMP-T ( I - 1 ) >/ { T ( I ) -T( I - 1 ) > 
301 CCNTINUE 
Rt TURN 
30? CCNTINUE 
3C3 CCNTINUE 
250 WRITE CG.2C5) ENTH 
STOP 

r no 


FUNCTION THCS(T) 

THIS FUNCTION COMDUTFS THT THERMAL CONDUCTIVITY OF SOLID THf RMKEEP 
AT 2 SCF, IT IS NOW A DUMMY FUNCTION WITH A CONSTANT VALUE 
IN KJ/M-SFC-r 
THCS=1 .0P00E-3 
RF TURF’ 

r nc 


FUNCTION HL(WL) 

This FUNCTION COMPUTES THC HEAT TRANSFER COEFFICIENT OF T HC 
HTATINC FLUID AS A FUNCTION Or FLOW RATE IN K J/SEC/M 2/K 
CRE=226P.576A 
RCCR IT=772A. 

CHTC=R. 91 3RE - A 
C V= . PF 03 
RE=WL«CRE 

IF(RE.LT.RECRIT ) GO TO AS 1 
HTCI=CHTC*RES$O.P 
Cf TO Ab3 
A5 1 hTCI = .C-tC9 

7 NU= 0 * lPl<=QX»(l.-n.B39/QX435.A/QX*<-?-2C7./0X**34A19./QX*«A) 
h TC I =2NU«UTC I 
AS 3 hL = H TC I 
RET URN 
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FUNCTION NT ( > ) 

THIS FUNCTION COMPUTES THE HEAT TRANSFER COE F F I C 1 1 NT OF 
VAPOR 1 7 I NG PQWf F T L U I D AS A FUNCTION OF QUALITY 
HT=X«.E01C74.3TR8A 
RETURN 
END 
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LIST OF SYMBOLS FOR MAIN PROGRAM 


I j L 
CL 

CS 

CT 

CV 

ET 

J6 

K5 

K7 


NC 

ND 

QC 

QL(I) 

QT 

TD 

T1 

T5 

T6 

WL 

WP 

ZL,ZP,ZPS, 

ZPV 

QCHG 

TIME ( I ) 


Integer loop counters 

Specific heat of transfer fluid, 
kJ/kg-°K 

Specific heat of subcooled power fluid, 
kJ/kg-°K 

Cycle time, seconds 

Specific heat of vaporized power fluid, 
kJ/kg-°K 

Total elapsed time, hours 

Indicator 1 = first time through 

2 = a new Thermkeep mass has 
been entered 

The number of entries in the table of 
solar collector output vs. time 

Indicator 1 = single fluid flow 

2 = transition between single 

and two-fluid flow 

3 = two fluid flow 

Number of daily cycles run 

Maximum number of daily cycles desired 

Rate of solar collector output, kW 

Rate of solar collector output at time 
I, kW 

Heat of vaporization of power fluid, 
kJ/kg 

End of the daily cycle, secs 
Start of the daily cycle, secs 
Beginning of net charge, secs 
End of net charge, secs 
Transfer fluid flow rate, kg/sec 
Power fluid flow rate, kg/sec 

Coefficients used to determine time step 

Heat rate storage, kJ 
Time which matches QL(I) 
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TLIN 

TPIN 

TSAT 

WMAX 

XOUT 

ASTOR 

DQCHG 

DQDIS 

DTIME 

DTIML, DUMP 
DTMAX jDTMIN 
QBURN 

QLOSS 

TAPPR 

TLOUT 

TPDES 

TPOUT 

PCDIFF 

QDEFCT 

QDISCH 

QEXTRA 

QPOWER 

QSOLAR 

TALWBL 


Transfer fluid inlet temperature, °K 

Power fluid inlet temperature, °K 

Power fluid saturation temperature, °K 

Maximum transfer fluid flow rate, 
kg/sec 

Power fluid outlet quality 
Amount of storage medium, kg 
Increment of storage input, kJ 
Increment of storage output, kJ 
Time step, secs 

Temporary time step calculations 

Maximum and minimum time steps, secs 

Auxiliary heat required to maintain the 
desired outlet temperature, kJ 

Solar collector heat unable to be used, 
kJ 

Transfer fluid outlet temperature 
approximation at beginning of two-fluid 
flow, K 

Transfer fluid outlet temperature, °K 

Desired power fluid outlet temperature, 
°K 

Power fluid outlet temperature, °K 

Percent different between storage heat 
in and heat out 

Heat required to raise the power fluid 
outlet temperature to the minimum 
allowable temperature, °K 

Storage output, kJ 

Heat from power fluid exiting above the 
desired outlet temperature, kJ 

Boiler requirement, kJ 

Solar collector output, kJ 

Minimum allowable outlet temperature, °K 
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LIST OF SYMBOLS FOR SUBROUTINE ELANAL 


E 

I,J,K 

T(I) 

Z(I) 

AI ,BI ,CI 
FI(I) 

A1,A2,A3, 
B1,B2 ,B3 

CC 

CH 

CL 

CS 

CT 

CV 

ET 

HL(WL) 

HT( X) 

J6 

K1-K8 
LI ,L2 
Ml 


The unit total energy content, kJ 
Integer loop counters 

The temperature of the Ith element, °K 

The solid fraction of element I 

Coefficients of unknowns in diagonals 
A, B, C used in subroutine TRIDIA 

The coefficient of the result of equation 
I used in subroutine TRIDIA 

Coefficients of unknowns used to calcu- 
late AI, BI, and Cl 

Product of the power fluid flow and 
specific heat 

Product of the transfer fluid flow and 
specific heat 

Transfer fluid specific heat, kJ/kg-°K 

Specific heat of subcooled power fluid, 
kJ/kg-°K 

Cycle time, secs. 

Specific heat of power fluid vapor, 
kJ/kg-°K 

Elapsed time, hours 

Function which returns the transfer fluid 
heat transfer coefficient as a function 
of flow rate 

Function which returns the power fluid 
heat transfer coefficient as a function 
of quality 

Indicator 1 = first time through program 
2 = restart cycle with a new 
storage mass 

Integer loop counters 

Integer loop counters 

Condition indicator: 1 = transfer fluid 

was calculated to increase going through 
an element 
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M 3 

ND 


Condition indicator: 1 = vaporization 

occurs completely in bottom element 

Number of daily cycle in progress 


NE , XNE 
QA 

QB 

QD 

QP 

QT 

TD 

TL(I) 

TP(I) 

T5 

T6 

WL 

WP 

XI(I) 

XP(I) 

AEL 

AEP 


Number of elements used in numerical 
analysis 

Rate of heat flux into the power fluid 
stream in an element above the previous 
time step's point of total vaporization, 
kW 

Rate of heat flux into the power fluid 
stream in an element be^ow the previous 
time step's point of incipient vaporiza- 
tion, kW 

Rate of heat flux into the power fluid 
stream in an element in the previous 
time step's transition region, kW 

Heat transferred to the power fluid in 
a whole element in the 2-phase region, 
kJ 

The heat of vaporization of the power 
fluid, kJ/kg 

The cycle ending point, secs 

Temperature of the transfer fluid at the 
lower boundary of element I, °K 

Temperature of the power fluid at the 
lower boundary of element I, °K 

Beginning of net charge, secs 

End of net charge, secs 

Transfer fluid flow rate, kg/sec 

Power fluid flow rate, kg/sec 

The solution to the Ith equation in 
TRIDIA, °K 

The quality of power fluid at the lower 
boundary of element I 

The total contact surface between the 
transfer fluid and the storage medium, 

m2 

The total contact surface area between 
tije power fluid and the storage medium, 
in 


D-36 



APL 

The total contact surface area between 
the power fluid and transfer fluid, m 2 

CCS 

Product of the power fluid flow rate 
and subcooled specific heat 

ccv 

Product of the power fluid flow rate 
and superheated vapor specific heat 

CT1 

Half-hour timer for display purposes, 
secs 

DQE 

The rate of energy change of an element 
during a time step 

EFF 

Heat exchange effectiveness 

ELM 

The storage medium mass in one numerical 
element, kg 

FTP 

The fraction of an element in which the 
power fluid is in 2-phase flow 

FXF 

The fraction of element NEF which con- 
tains less than quality 1.0 fluid 

FXI 

The subcooled fraction of the element NEI 

NEF 

The element in which the power fluid 
reaches quality 1.0 

NEI 

The element in which the power fluid be- 
gins to vaporize 

QEC 

The heat rate from storage medium into 
the power fluid stream 

QHC 

The heat rate from the transfer fluid 
into the power fluid stream 

QHE 

The heat rate from the transfer fluid 
into the storage medium 

QTP 

Heat transferred to the power fluid in 
the two-phase region of an element, kJ 

RES 

Resistance to heat transfer into the 
power fluid 

RXF 

The fraction of element NEF in which the 
power fluid has greater than 1.0 quality 

RXI 

The fraction of element NEI in which the 
power fluid has greater than 0 quality 

SUM 

Anticipated heat to be transferred to 
the superheated power fluid, kJ 

TEL(I) 

An element temperature from the previous 
time step 
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TLP 


TLI 

VAP 

XII , XI2 ,XI 3 

XLE 

XNE 

XNT 

XPA 

XPB 

XP1 ,XP2 

AELE 

AEPE 

APPE 

APLE 

CPST 

DEPF 

DRXP 

DSTK 

ENEL 

ENER(I) 


The temperature of the transfer fluid at 
point RXF, °K 

The temperature of the transfer fluid at 
point PXI , °K 

The vaporization subroutine 

Temporary points used in an iterative 
search for the point at which vaporiza- 
tion ends 

The length of a tube in an element 

The number of elements used for analysis 

The number of tubes in each heat ex- 
changer 

The calculated power fluid quality leav- 
ing an element which contained the end 
of transition the previous time step 

The calculated power fluid quality leav- 
ing an element which completely contained 
transition the previous time step 

Temporary qualities used in an iterative 
search for the point at which vaporiza- 
tion ends 


Contact surface area between storage med- 
ium and transfer fluid in one element, 

2 

nr 


Contact surface area between storage med- 
ium and power fluid in one element, nr 


Inside surface area of the power fluid 
tube in one element , m^ 


Contact surface area between transfer 
fluid and power fluid in one element , m^ 

The specific heat of steel, kJ/kg-°K 

The effective diameter of dual heat ex- 
changer tube 

The calculated change in RXP 

The diameter of solid Thermkeep buildup 
on a tube, m 

The energy content of one element re- 
turned by subroutine TKPROP, kJ 

The energy of the Ith element, kJ 
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FRAC 

FSUP 

HRTS 

HRT6 

HTCL 

HTCS 

HTCV 

QANT 

RHOM 

STME 

STMM 

THCS 

TLLN 

TPIN 

TSAT 

UAEC 

UAHC 

UAHE 

VSTK 

XNTU 

XOUT 

YNTU 

ZNTU 


The solid fraction of an element re- 
turned by subroutine TKPROP 

The fraction of an element which con- 
tains superheated power fluid 

Elapsed time at time T5, hrs 

Elapsed time at time T 6, hrs 

Transfer fluid heat transfer coefficient, 
kW/m2-°K 

Heat transfer coefficient of subcooled 
power fluid, kW/m 2 -°K 

Heat transfer coefficient of power fluid 
vapor, kW/m 2 -°K 

Heat into the power fluid anticipated 
for section FXF, kJ 

Storage medium density, kg/nr 

Steel mass in an end element, kg 

Steel mass in a middle element, kg 

Thermal conductivity of solid storage 
material, kJ/m-°K 

Inlet temperature of transfer fluid, °K 

Inlet temperature of power fluid, °K 

Power fluid saturation temperature, °K 

Conduction coefficient of storage med- 
ium to power fluid, kW/°K 

Conduction coefficient of transfer fluid 
to power fluid, kW/°K 

Conduction coefficient of transfer fluid 
to storage medium, kW/K 

The volume of solid storage material in 
one element, nm 

Number of heat transfer units in an ele- 
ment 

The power fluid outlet quality 

Number of heat transfer units in the two- 
phase section of an element 

Number of heat transfer units in the 
superheated section of an element. 
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ASTOR 

CHGTM 

CRKTL 

DTIME 

ETSAT 

QCRXF 

QCRXI 

QPART 

RESTP 

RSOLL(I) 

TINIT 

TLOUT 

TPOUT 

UAECF 

UAECI 

UAEC2 

UAHCF 
UAH Cl 
UAHC2 

UAHEF 


Amount of storage material, kg 

Elapsed time from beginning of net 
charge, hrs. 

Intermediate variable = 

2tt x XLE x XNT used in calculating 
solid layer resistance 

Length of time step, secs 

Intermediate variable used to calcu- 
late YNTU 

Heat into the power fluid in section 
RXF, kJ 

Heat into the power fluid in section 
RXI , kJ 

Heat into the power fluid in the two- 
phase section of an element 

Resistance to heat transfer into the 
power fluid in the two-phase region 

The solid layer resistance in element 
I 

Unit starting temperature, °K 

Transfer fluid outlet temperature, °K 

Power fluid outlet temperature, °K 

Conduction coefficient of storage medium 
to power fluid in section FXF, kW/°F 

Conduction coefficient of storage med- 
ium to power fluid in section RXI, kW/°K 

Conduction coefficient of storage med- 
ium to power fluid in a purely transi- 
tion element, kW/°K 

Conduction coefficient of transfer fluid 
to power fluid in section FXF, kW/°K 

Conduction coefficient of transfer fluid 
to power fluid in section RXI, kW/°K 

Conduction coefficient of transfer fluid 
to power fluid in a purely transition 
element, kW/°K 

Conduction coefficient of transfer fluid 
to storage medium in section FXF, kW/°K 
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UAHEI 

UAHE2 

XLEPF 

XPXI3 

CUAHCS 

CUAHCT 

CUAHCV 

NEFTMP 

NEITMP 

OLDQEC(I) 

QSUPHT 

RESSUP 

TTRIAL 

XNTUSP 


Conduction coefficient of transfer fluid 
to storage medium in section RXI , kW/°K 

Conduction coefficient of transfer fluid 
to storage medium in a purely transition 
element 

Tube length, m 

Intermediate variable used in an itera- 
tive search for the vaporization ending 
point 

Conduction coefficient of transfer fluid 
to power fluid in partial elements for 
subcooled, transition, and vapor phases, 
kW/°K 

Element where end of transition is 
temporarily located 

Element where the beginning of transition 
is temporarily located 

Heat transferred from storage to power 
fluid in the Ith element during the 
previous time step 

Heat into the power fluid in the vapor 
region of an element 

Resistance to heat transfer into the 
power fluid in the superheated region 

A calculated power fluid temperature at 
the beginning of two-phase flow, °K 

Number of heat transfer units in the 
superheated fraction of an element 
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APPENDI X E 
SINGLE HEAT EXCHANGER COMPUTER CODE 
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C MA IN PROGRAM TESPR 

COMMON CPF .7U.TL.XNE » PM AX • PDM AX .tofPUMP.tofPt AK . TF ( E J ) ,CT .AST OR ,EN 3 »WC 
COL L 

D I HE NS 1 OF! GOUT (25 ) . T IMC 25 ) 

30 G FERMAT ( I 3.6F It . 1 ) 

301 FORMAT! 1 9X • * 7ABLC Of SOLAR COLLECTOR OUTPUT VS. T I HF •) 

3Ci.‘ F t RMA T ( 24 X • * PATE IN KW *,9X.» TIME IN SECONGS *) 

303 FCRMA7(F 1G.4 .F10. 1 ) 

30 A FORMAT <2F 10.4 ,2F 10. 1 .F1G.0.F10.2.F1C.3) 

30 E F t RMAT ( 5X • * FLUID SPfCIFIC HEAT. KJ/KG/K •.rjC.A) 

306 F t PMA T ( 5X » • TF-ERMKCEF SPFC HEAT. KJ/KG/K *.F10.4) 

307 F C P'aAT < 5X . * TURBINE INLET TEMP.* K '.F1C.1 ) 

30 E F t RMAT ( OX . * TURBINE OUTLET TEMP., K *«F10.1) 

3C9 F t RMAT (EX,* NO. OF E l E ME NT S FOR ANALYSIS •.FIO.C) 

310 F C RMA T ! 5 X * * AFT. OF STORAGE MAT L. , KG *.F1G.2) 

311 F C RMA T ( 5X , * MAX. PUMP FLOto , KG/SEC. *.F10.3) 

31? F C RMA T f 5 X « • MAX. TIME INCREMENT. SFC. *,F10.2) 

313 F 0 RMAT ! 5X » • BI1LER DEMAND. KW *.riG.3) 

31A FERMAT ( 1H1 .EX. * MAIN PROGRAM DATA*.//) 

3 1 5 FCRMAT(25X.F JC.4.9X.F10.1 ) 

316 FCRMATCEX,* NO. OF DAILY CYCLES *,7X.13.///> 

317 P t PMAT ( 1H1 » 5 X » * HEAT OF FUSION THERMAL STORAGE ANALYSIS*.///) 

320 FORMAT (///,?0X»* CYClf DATA*.//) 

321 F OPMAT ( 1 OX » » ELAPSED TIME. HR «.F10.1) 

322 FCRMATIJOX.* CYCLE TIME, SEC *,F1G.2) 

323 F E RMAT ( 1 GX . * LOST COLLECTOR HEAT, KJ *,E15.7) 

32 A FERMAT! 10X,* CCLLECTCR OUTPUT. KJ *»Elt.7) 

32 5 F C RMAT ( ] OX , • STORAGE OUTPUT. KJ *.E15.7) 

326 F t RMAT ! 1 OX . • STORAGE INPUT. KJ *.EJ5.7) 

327 FLPMATI1CX.* AUXILIARY INPUT, K « *«E15.7) 

32E FERMAT! 10X.* ECILEP REGUIPEMENT. KJ *.E15.7> 

329 F C RMAT C 1 OX » * FEAK STORAGE FLOW, KG/S *.FIO.A) 

330 F C RMAT ( 1 OX » • tEAK PRESS DROP. KCF/SM *.E1C.4) 

334 F C RMA T t 1 OX , * PEAK PRESS DROP. PSI *,E10.4) 

331 F C RMAT f 1 0 X • • TEAK PUMP PChER. Kto * .E 10 .4 ) 

332 FCRMATI10X,* TOTAL PUMP WORK. KW-HR *.E10.4) 

333 FORMAT ( 3FJ0.2 »F 10 .3. 13) 

33 E PEP MAT! 1 0 X * * NUHEFR CF TIME INC RE ME NTS * . 1 1 0 . /// ) 

33C F C RMAT ( S X , • MIN. T I MF INCREMENT, SEC. *.F1C.L) 

337 F t RMAT (EX.* IDLE TIME INCREMENT. SFC. *.FIP.2) 

336 F C PMAT I 1 OX* * DEFICIT REQUIREMENT, KJ * »E IC .4 , 1 2> • • CYCLF TIME 
tFLCW Tf ( 1 ) TFINE + l) 7F(1)-F TF IKE •* 1 >-F * ) 

340 FCPMAKF 10. 0) 

RL AD! 5. 300 ) K E .7 1 .7 2 .73 ,T4 .T5 .7 6 
PE AD! 5. 304 ) CFF.CTK.TU.TL. XNE.ASTOR.HMAX 
Rf AD ! 5 » 333 ) LTMAX.DTF 1N.DT1DLE. GE »ND 
Jf =1 

36C CONTINUE 

IF (J6.EG.2 ) GL TO 37f 
to R ITE ( 6 .3 17 ) 

HR ITE (6.30 1 ) 

HR ITE16.3C2 ) 

DC 1 1-1. KE 

RE AD! 5 . 303 ) GOUT ( I ) , T 1 M ( 1 ) 

1 hR ITE(6.315 ) GOUT ( I ) .TIM! I ) 

370 HR ITE (6.314 ) 
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KP I TE < 6 • 305 1 CPF 
WP 1 TE ( 6 • 306 ) C7K 
W P 1 TE C 6 * 3C 7 1 TU 
KP I TE t 6 * 308 1 7L 
WP 1 TE < 6 .309 ) XNE 
*P ITE (6 .310 > ASTOR 
KP1TE(6.311) PP'AX 
NP 17E(6»312) CTMAX 
KP1TEI6.336) DTM1N 
WP I TE ( 6 .337 ) CT IDLE 
KP1TEC6.313) GB 
W P ITE (6.3161 ND 
tC 1 = 1 

DHMAX=CPF*( T tl-TL 1 
WM IN=GB/DHMA> 

ET=0. 

OC EF CT = 0 • 

GBClL=0P*CT2-*7e-7 3) 
hr P EAK =0 . 

PM AX=0. 

WPUMP=0. 

PL MAX-0 . 

ND 7=0 
NC = 1 

C INITIAL J2A7 ION DF DAILY CYCLE 

C K7J 1=N0 FLO*. 2=CHAPG1NG. 3=D 1 SCHARG I NG 

AC C7=DTMAX 

£7=ET4DTMAX /3( 00 . 

X 7 = 1 

QLCSS=0 . 

GSCLAR=GB«DTMAX 
QD I 5CH- 0 . 

GCHG=0. 

QE URN = 0 • 

W F =0 • 

DT 1ME =D7MAX 
QC =GB 
L= 0 

50 L=L+1 
NL7=NDT-» 1 

CALL EL ANAL ( V.I .T1.T2.73.T4.T5.T6.K7.DTIME.K1.JC) 

IF(L.NE.IC) GC TO 51 
7 E ML= 1 • 6 *T F ( 1 >-4 59.6 7 
7 E MU= 1 .6 $T F CNt 41 1-459.67 
L = 0 

51 Z- .06SAST0R*C7K/(XNE *CPF 1 
IF (K7.EG.1) GE TO 60 

IF (K7.EG.3 ) GC TO 55 

GL 0SS=GL0S5> + <f,C-QB“toM-CPF$ (TU-T FC 1 1 1KDTIME 
Gi OLAR = GSOLAR-*GC6DT I PE 
QCHG=GCHG4KF*CPF* (7U-TF C 1 1 )*DT1 ME 
GC TO 7C 
55 NE =XNE 

QECRN=oeURN+*r#CPF=t(7U-TF<NE41 1 )$DT I ME 

IF CTFCNE + l 1.LT.5B0. ) GDEF C T =GDt FCT 4 *F*CPF « ( 560 . -TF ( NE 4 1 ) } *D T I ME 
GD 1SCH=GD1SCH4 |i;F $CPF «(TFINE41)-TL 16DTIMF 
CC» G50LAR = GSCLAP4CC*D7 IPE 
70 DC 2 1=2. K5 
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IF (CT.GE.TIMd )) GO 1C 2 

GC = QOUT € 1-1 )-* ( QDU1 ( 1 l'GOUT (1-1)>*CCT-T1M(1-1))/(T1M(I»-TIM(1-J>> 
GC TO 60 
2 CONTINUE 

C INADEQUATE COLLECTOR OUTPUT AT END CE DAY. EVEMM 

eo lF(CT.CE.Tl) GO TO 9C 
KCCLL-GC/DHMAX 
M F =WMIN-MC0LL 
K7=3 

D 1 1ME=2/KF 

IF (DT1NE .LT.DTKIN ) DT1ME=DTM1N 
IF (DT1ME .LE .DTKAX 1 GO TO 61 
D 1 IME=DTKAX 
FI CT=CT+DT1ME 

IF (CT.LE.T1 > CO TO F2 
DT IME=DTIME-(CT-T 1 J 
Cl =T I 

F2 El =FT-*DTlME/3tC0. 

GO 10 50 

C NO HEATING AVAILABLE FROM COLLECT ORE 

90 IF (CT.GL.T2) CL TO 1 TO 
*F =KMIN 

KC 0LL=0. 

K7 =3 

DT I ME = Z/WF 

IF (DT IKE .LT.DTMIN ) DT1ME=CTwIN 
IF (DTIME .LE .DTMAX ) GO TO 91 
DT I ME = D T MAX 

91 C T =CT ■•DTI ME 

IF <CT .LE .12 ) CO TO 92 
DT ]ME=DTIME-(CT-T2) 

CT=T2 

92 FT =ET+DT 1ME/3F 00 . 

6C TO 50 

C OVERNIGHT IDLf 

IOC IF(CT.GE.T3) CC TO I 1C 
*F =0 
N7 = l 

DT 1KE-DT1DL E 
C T =CT T DT I ME 

IF (CT.LE .T3> CO TO 1(2 
DT IME=DTIME- (CT-T3) 

CT =13 

102 ET=ET-»DT 1ME/3C00. 

GC TO 5C 

C NC HEATING AVAILABLE FROM COLLECTORS AFTER OVERMCHT ILLE 

11C IF(CT.GF.TA) CO TO 120 
VF =*M1N 
K7 = 3 

DT IME=Z/KF 

IF (DTIME .LT.DTMIN ) D 11 ME =D T MIN 
IF (DT IKE .LE .DTMAX 1 GC TO 111 
DT IME =DT MAX 

111 C T =CT +DT 1 ME 

IF (CT.LE. TA ) CC. TO 1 1? 

D 1 IME=DT1ME-(CT-TA ) 

CT =T A 

112 ET =ET+DTIME/3CC0- 
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Gt TO 60 

C INADEQUATE COLLECTOR OUTPUT AT E'f GINNING OF PAY, FORKING 

12C IF (CT.GE .75) CO TO 1 ?C 
AC C'LL = QC/DHM A> 

AF = AM I N” A COLL 
K7 = 3 

DT JME = Z/Wf 

IF (DT1ME .LT.DTNIN > DT1ME=DTWIN 
IF (DT1ME .LE .L'TMAX ) GC TO 121 
DT 1ME=DTMAX 
121 C T -CT ■♦DTI ME 

IF CCT.LE.T5 ) GO TO 122 
DT 1ME=DTIME -CCT-T5 ) 

CT =T6 

12? ET=ETTD7IME /3< 00 . 

GC TO 50 

C EXCESS SOLAR HEATING AVAILABLE FOR STEKAGL 

13C IF C CT. GE . TG ) GO TO 1 AC 
IF C CT .GT.T5) CO TO 1 31 
AF =0 
K7 = 1 

DT IME=DTMAX 
GC TO 132 

131 AF = CQC-GB I/CCrF *C TU- TF C 1 1 1 ) 

ACCLL=#F 4AM1N 

IF (ACCLL.GT.AMAX) WF=WMAX-AVIn 
IF CVCOLL .GT .V.FAX ) WCCLL=AMAX 
X 7 =2 

DT 1ME=Z/AF 

IF (PT1ME.L T.0TM1N ) D T I ME =D T M IN 
IF COT JME .LE .DTMAX 1 GC TO 132 
DT IME=DTMAX 

132 C T =CT -FDT 1 ME 
IFCCT.LE.T6) CO TO 133 
DT IME-DT1ME-CCT-TC) 

CT = T6 

133 ET=ET+DTIME/3f 00. 

GC TO 50 

C CUMULATIVE INFORMATION AND MAIN F-Pf GRAM END 

1 A 0 A R ITEC6* 32 0 ) 

AR 1TEC6.321 ) t 7 
AR 1 TE (6.323 ) GLOSS 
ARITEC6.32A) GSDL AR 
AR17EC6.325) GDISCH 
MR 1 TE C6.326) GCHG 
AR JTEC6.327) GE UR N 
AR ITE (6. 338 ) GDEFCT 
AR ITE C6. 32 6 J GBDIL 
AR ITE (6.329 1 V PE A K 
ARITEC6.33C) F DMA X 
PS IMAX=PDMAX*.001 422 
ARITEC6.3341 FSIVAX 
AR ITEC6.331 ) F M AX 
A fi 1 TE C 6 .332 I APUVF 
ARITEC6.335) NDT 
NC =NC+1 

QC IFF 1=0.5*EN? 

GD IFF?=0CHG-GI 1SCH 
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IF <GDIFF?.GE .CDIFF 1 ) CD TD 99 
lFCNC.GT.ND) CD TO 99 
Q L CSS =0 » 

G £ OL A R = G • 

Q t I SCH=G • 

OC HG=0. 

GBURN=0 . 

ODEFC T = 0 • 

MrXAIC = 0. 

PL MAX=0 . 

PM AX=0. 

WPLMP=0. 

ND T =0 
GC TO 40 

99 R E ADC 0*340) A£ TOR 

IF (AS10K.EQ.1.) CO TT 360 
JC = 2 

GO TO 360 
360 STOP 
END 


.. f : 
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SLBROUT 1NE ELANALCwF ♦ T 1 .72 • T3 . T 4 . T E .T6 . K7 . DI IKE » M » J6 ) 

C ThlS SUBROUT1FF PERFTPMS THE HEAT TRANSFER COMPUTATIONS OK THE 

C ELEMENTS 

DIMENSION EN£ GC5 1 ).2CM).T(51).RMAC51) 

CCMMON CPF .TU.TL .XNE .FMAX.PDMAX .WPUKP.WPt AK. TF CL 1 ).CT,AS7DR.EK3.».C 
60LL 

CLMMON/PROPTK/ELM .CPST.STMM.STME 
IF (J6.E0.21 CC TO 3BC 
GC 10(1. 2). K] 

C FIRST TIME THROUGH - INITIALIZE AND READ DATA 

1 NE=XNE 
CT 1=0. 

RE ADCS. 101 ) DC.DI .XL T . XNT . T INS . COO . D ID 

RE AD C E> • 1 02 1 TEND. TSHf LL .DHEL 

RE ADC 5. 103) RHOST .RHOF .RHOM.RHO INS 

RF AD(b. 104 ) THCI NS.THCST .1 HCF.PRF »V ISCF * CPST 

R E AD 15.105) TAMB.TSHR.T INI T 

RE ADC 5. 700 PKATS »PF Ef A *PF SSW.PFSTA.PFSTX.Pf SB A .F f SEW 
RE ADC 5.701 ) GA T OPS . G A S I DE . GABOTS .PMVW.PM VE .PF VV .RFV* 

RE ADC 5. 702) CPMTW.PFTX.PFTA.PTK. PINS. ASPECT 
360 CONTINUE 

IF CK1.EG.2) GC TO 2 
C GEOMETRICAL CALCULATIONS 

V7UBES=.76544>LT4XNT4C0442 

VEL7K=AST0R/RHCM 

VTOT=VOLTX+VTLEES 

t ATTEMPT SPECIFIED ASFFCT RATIO 

JZ =1 

DSHCLL= CVT0T41.27324/ASPEC T )44( 1./3.) 

HE HELL=DSHELL4 ASPECT 

IF CHSHELL.LE .XLT ) GO TO 3 

HSHELL=XL7 

DSHELL=2.4S0KTCV70T/C3.141 64HSHELL ) ) 

JE =2 

2 DSHR=DSHELL+2.4T 1 NS+ 2 .4TSHELL 

C WEICHT CALCULATIONS 

HGTT0T=HSHELL+2.471Nf*2.4T END 
KTENDS= I .S7064TEND4RFCST4D SHELL 442 
Ml SHEL=3. 14 1 64DSHELL47SHELL4HSHCLL4RHDST 
*TVESS-*TSHEL*X TENDS 

KT TUBE =.7654* (DD4«-( 1«*2) #RH0STtXLl*XN1 

MT INS=RH01NS4< . 7 O 54 * C P SHR4 *2 - C D SHEL L* 2 . 4 T E ND ) 4 -4 i ) 4 HSHE LL + 1 . 6 70 B 4T I 
4NE 4DSHR442 ) 

C COST CALCULATIONS USING SIX-TENTHS FACTOR 

ATUEiES = 3.14 1 C 4 XL T 4XK T < DC 
PMTk=CPMTX/DC 
CMTUBE = MTTUBE 4PM TW 

CC TUBE = CMTUBF *23.61 6 24PFTX 4 C WT T LP E / 4 7 .2326 ) 4 4 . 6* i . Of 7 4PF T A 4 C A TUBF S 
4/2 .1343 ) 44 • 6 
C0ST0R=AST0R4PTK 
CC INS = M TINS4P1NS 
ASIDES=3.14164DSHR4HCTT0T 
XTS1DE = AS1DES4GAS IDE 4RH0ST 
CMSIDE=WTSIDL4PMA TS 

CCS IDE = CMS IDE *23.92 4FFSSM4 CXTSI0E/47.F4) 44.6*2.1 41 FLSA4C AS1CE5/5. ) 
444 .6 
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AeOTS = .7e5A$E/SHR**2. 
hTECTS = Aet)7S#CABD7£ A f HOST 
CPPCTS=»*TeCTS<PMA7S 

C0B07S = C MOOTS 46. I7=tFr£BW#(W7B07*-/lC.34)4*.e-«.3C5 C5*F FSPA * f AE 07S/* C 
91 19)#«.6 
A7CPS=AE0TS 

WT 7CPS = AT0PS*GAT0PS9RH0ST 
CPTOPS=bTTOPSSPMATS 

CC70PS = CM7 0Pi>-*3.f t*Pf STW9( W 770P S/7 . 32 ) *9 .0 4 . 30 1 959 PF S7 A* C A 7 CPS / . t 1 
C19 >99.6 

CFi>HRD = CMSIDL-»CMB07S-»CM70PS 
CC. SHRD = C0£IPt" , C0BO7S-*C07DP5 
CP VESS=WTSHEL9PMVP4t»TrNDS9P¥VE 

CC VESS = CMVESS-»90 .59 PE VK9(b TVESS/ 16 1 . } 99 . C+ .26 A 9C T V V 9 i V7C7 / . BCE ) 99 . 
*6 

CP 70T=CMTUBE4CFSHRD+CPVCSS 4C0INS4C0ST0R 
CC707=C07UBC-»C05HRD-» CTVE SS 4 CO IN S4 C0ST0R 
CC FTES=C0T07 -CP70 7 
Vt 7 SHRD = W TSlDE4PTT0PS4p7bCT S 
b7 707=b7VESS+V77UPE4V7INS-» AS70P *V 7 SHRD 
( RESISTANCE CALCULATIONS - LA7FRAL 

FI AREA=0.785A9>N7«DI992 
RE CRIT = 2100 .9( 1.4 J?. * <D 1/BHEL >990.5 ) 

CRE=DI/EFLAREA9V1 SCF ) 

CF 7C=0.0229 7HCF9 (PRF 99 .4 )/Dl 
AI UBE1=3.141C9XNT9DI9XLT/XNE 
CRKTL=6.2e39XL79XNT/XNE 

RINSL=ALOG (DSFR/DSHE LL ) 9XN E / ( 6 . Z 6 3* 7 HC I NS9HSHE LL } 

C RESISTANCE CALCUL A7 I CNS - AXIAL 

R7UBEA=XLT/( .785 A9THCST9XNT9XNC 9 < D099 2-D I 992 ) ) 

RS HEL A=H£ HE LL/ ( 3 • 1A 1 C9XNE9 7HCS7 $7 SHELLED SHELL) 

CPKMA=HSHELL/( .785A9XNE9EDSHELL 992-XNT 9D0992 ) ) 

C RFSIS7ANCE CALCUL AT 1 f NS - END 

R7UBEE=7 INS/C .785 A*7» CST9 ( D00992-D I 0992 ) I 
R 1 NSE =7 1 NS/ ( .7e5A9THClNS9{ DSHR992-XNT9D0992 ) ) 

RAXE=1 ./( 1 ./R7UBEE-* 1 ./RINSE ) 

( MASS CALCULA710NS 

ELP-AS70R/XNE 

£T PM = 3. 1416*RF0ST9( DSHELL9 T SHELL 9HS HEL L4 XL T 9XNT 9 ( E 0**2 *D I 992 )/ A . ) / 
1XNE 

ST ME =STPM4 .7t549PH0S T 9TEND«(DSHELL992-XNT9Dl992 ) 

1C1 F C RMA TI2Fl0.e»F10.2»F10.0»3F10.5) 

102 F 0 RMA T ( 3F 10.6) 

103 FC RMATE 4F 1 0.2 ) 

1 OA FE RPAT C3F 10. WF 1 0 . 5 » f 1 0 . 9 » F 10.5 1 
105 FC RMATE 3F 10 .2 ) 

IOC FORMATE///. 12X»* SUBROUTINE ELANAL DATA*./) 

107 F E RMAT E SX » • SHROUD DIAMETER. M. »»F10.4) 

10e F E OP A 7 ( EX » • VESSTL DIAMETER. M. •.F 10 .A) 

109 F CRMAT ( 5X . • TUBE O.D.. M. ••F1C.6) 

110 FORMAT E 5X » • TUBE I.D.. M. ••FI0.6) 

111 F L RMAT ( 5X • * VESSEL HC 1 GHT . M. »,riO.A) 

112 FCRPATE5X.* TUBE LENGTH. P. »,F10.2) 

113 F C RMA T ( 6X . • NC. OF TUBES ••FIO.C) 

1 1 A F0RMATE5X.* INSULATION THICKNESS. “ . ».F10.5> 

115 F D RMA T ( 5X . * VESSEL P ALL THICKNESS. M.»,F10.6) 
lit F C RMA 7 ( 5X » • VESSEL E HC THICKNESS. P. *.F10.C) 

117 F 0 RMA T ( EiX » • THERMKEEP DENS., K G/C U .P . • . F 1 0 . 2 > 
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F C R M A T C 5 X * 
FI RMAT (EX , 
FC RMAT ( EX » 
Ft RMAT (EX, 
F C RMA T ( £ X * 
FORMAT (EX, 
FCRMAT ( 5X • 
FORMAT C5X. 
FORMAT (EX. 
F CRMAT ( EX » 
FC RMAT(EX, 
FC RMAT! 10X 
FCRMAT (EX . 
FCRMAT (EX. 
FCRMATC 5X» 
FCRMATC6X. 
FC RMAT (EX. 
FCRMAT ( 6X . 
FCRMAT (EX. 
FCRMAT(EX. 
FORMAT! 1 HI 


FORMAT! 
FC RMAT ( 
FORMAT! 
FORMAT! 
F C RMA T ( 
FORMAT! 
Ft RMAT ( 
FORMAT ( 
FORMAT! 
FORMAT ! 
FORMAT! 
FCRMAT ( 
FORMAT ( 
FORMAT! 


STEEL DENSITY. KG/C l .M • 
FLUID DETflTY, KG/CUM. 
TF. 0 OMC . 1 NSUL.. KW/M-K 
TF. C OND . STLEL. KW/M-K 
TO. COND. FLUID. KW/M-K 
FLUID PRJNDTL NUMEER 
FLUID VISCOSITY. KC/M-SCC 
ATE'lENT TEMP.. K 
SOROUD TFMF..K 
IK 1 T 1 AL T C MP . • K 
STEEL SP. HT.. KJ/KG-K 
' ASPFCT RATIO NOT MET FOR 
INSUL. DENSITY. KG/CU.M. 
SYSTEM HEIGHT. M. 

SYSTFM *C IGHT * KG. 

VC SSE L WC 1 GOT . KG. 
hi IGHT 01 TUBE S * KG . 
INSULATION WEIGHT, KG. 
CLTLrT 0 • 0 . , M. 

OUT L E T 1 .D. . M . 

' EEGINN 1 NG OF DAILY CYCLE 


INTER ME D I A T E 
INTER ME DIATF 
INTtRMED 1ATC 
INTERMEDIATE 
I NT! P MED I A T t 


KG/CL. M. *.FK.2) 

KG/CL. M. *,F10.2) 

L . . K W/M-K • .T 10.E ) 

L • KW/M-K • .F 10. fc I 
D. KV/M-K ••FlO.fc) 

NUMEER *.F 1 G • E ) 

Y. KG/M-SE C* » F 1 G . S I 
K *.Fl0.2> 

• »F 1 C .0 ) 

K ••FIG.?) 

KJ/KG-K • *F 1 0 .t) ) 

NOT MET FOR TUEE LENCTl OF • 
. KG/CU.M. *»F1C.2) 

M. • » F 10.4) 

KG. • .F 1C . 1 ) 

KG. * »F 10.1 ) 

S. KG. * » F 10.1) 

GHT, KG. * »F 10.1 I 

. • .F 1C.S ) 

. • .1 10.5 ) 

DAILY CYCLE*.//) 

AT CT=T1 = *,F 10.1 .* SECS*.//) 
AT CT=T2= * .F 10. 1 , • SECS*.//) 
AT CT=T3=*.F10.1 .* SICS*.//) 
AT CT=T4=».F10.I .* SECS*.//) 
AT CT=T5= • .F 10. 1 . * SECT*.//) 


.F ir.2 ) 


END-CF-CYCLT SYS7C« MAR AT CT=T 6= * . F I 0 . 1 . * 


SE CS* .//) 
T(l)-F *.//) 


(1) T ( I ) 2(1) TF ( J ) T(l)-F *.//) 

X.I3.1X.F9.2 »3X»F9.2»3X.F9.2t3X»F9.2 ) 

X.I3.1X.F9.2.3X.F9.2) 

X,* T(TAL SYSTEM ENERGY, KJ • ,E 15*7 * /// } 

X.I3.25X ,F9 . 2 ) 

HI.* STARTING TOTAL ENERGY THIS CYCLE. KJ. *.E10.f.///) 
OX.F 1C.2 . 10X ,r 10. 1 , 10X.F10.5.10X.F 10 .2. I C> .F 1 0.5. I OX.F 1 0.2 


£ I 

FC RMAT ( T15 . *ET * .T 25 . ‘CHARGE TIME* »T45»*ST0R/ GO F LC W * .T 65 . • 7 -6fc OUT 
£• ,Te5,*C0LLECT0R F L CF * . T 1 0 5 » *RE TURN COME 1 NL D TEMP.*,/) 

Ft RMAT ! /,5X. I 1.5X.F lC.?«5X,F10.2,SX*E12.4»£X,t 12.4.fX«F12.4.EX»I2) 
FLRM,AT(7F10 .4 ) 

FC RMAT! 3F J 0.7.4F I 0. 4 I 
F t RMAT (OF J0.4 ) 

FCRMAT ( IH1 .3CX .* COST CALCULATICN INPUT DATA*.//) 


7E1 F CRMAT ( 30X , * SHROUD RATER] AL COST, t/F6. 

7 £2 F 0 RMA T ( 3 OX . * SHROUD SIDE FABRICATION COST. $/SG.f . 
7E3 F t RMA 7 ( 30X . * SHROUD SIDE FABRICATION COST. 1/K. 

7E4 F C RMA 7 ( 30X » * SOROUD TCP FABRICATION COST. S/SC.N. 

7SE F 0 RMA T ( 30X . • SHROUD TTP FABRICATION COST. S/KG. 

7E6 F C RMA T ( 30X . * SHROUD TCTTOM FABRICATION COST. */SG.M. 
7E7 F ORMA T ( 30X . * SHROUD BOTTOM FABRICATION COST. T/KC. 
7it F C RMA T ( 30X . * SHROUD TCP THICKNESS. M. 

7£R F CRMAT ( 30X. • SHROUD SIDE THICKNESS. M. 

760 F C RMA T ( 30X . • SHROUD FCTTOM THICKNESS. M. 

761 FCRMAT(30X.* VESSEL FALL MATERIAL COST. 1/KG. 

762 F C RMA T ! 3 0 X . • VESSEL END MATERIAL COST. S/KG. 

763 FCRMAT(30X«* VESSEL FABRICATION COST, S/CL.M. 

764 FORMAT (30X.* VESSEL FABRICATION COST. S/kG. 

765 F C PMA T ( 30X • * TUBE MATERIAL COST. $/K. 

76 f F ORMA T ( 30X • * TUBF FAORICAT I ON CCST. 1/KG. 


*,F 10.4) 
• • F 10.4) 
*.F 10.4 ) 
*.r ic .4) 
*.F I0.4I 
*,F 10.4 ) 
• » F 10.4) 
*.r io.7) 
*.F 10.7) 
*.F 10.7) 
*.F 10.4) 
• • F 10.41 
*.F IC.4 ) 
*,F 10.4) 
*,F 10.4 ) 
».F 10.4 J 


E-10 



767 F C PMA T 4 30 X • * 
766 FORMAT! 30X. • 
769 F I’RMA T 4 3 0 X. * 
77C FOPMAT45X.* 
771 F C RMAT 4 / • 30X 
77 ? F C RMA T C 30> * • 

773 F C RMA T ( 30 X • * 

774 FORMAT ( 30X > • 
77E FCRMATC30X.* 

776 F C RMA T 4 3GX • • 

777 F 0 RMAT 4 3CX • • 
776 FCRMAT(30X.» 
779 F C RMA T ( 30X * • 
7ec F 0 RMAT I 30X • • 
761 FORMA T ( 30X • • 
76? F 0 RMA T C 3 OX • * 
763 F C RMA T 4 30X * • 
76 A F ORMAT ( 30X • • 

765 FORMAT C 30X • • 

766 F OR’MA T ( 30 X • • 
7e7 FORMAT 4 30X.* 
76B FORMAT C30X.* 
769 FORMAT!/. 30X 
79C F ORMA T ( 30X • • 
791 F CRMAT ! f>X . * 


TUBC F A6R1CAT 1 ON COST. J/SG.M. 
STORAGE MATERIAL COST. T/KG. 
1NSULA71CN COST. WKG. 

ASPECT RATIO. HGT/LMA. *»F10.4) 

•• COST CALCULATION OUTPUT DATA*./) 
TUBE V A TT R 1 AL COST. * 

TOTAL TUPE COST. 1 
STORAGE MATERIAL COST. $ 

INSULAT ION COST, t 
SHROUO SIDE MATERIAL COST. S 
SHROUD SIDE TOTAL COST. 1 
SHROUD PCTTOM MATERIAL COST. * 
SHROUD EOT TOM TOTAL COST. S 
EHROUD TTP MATERIAL COST. * 

SHROUD TOP TOTAL COST* * 

SHROUD TOTAL MATERIAL COST, i 
SHROUD TOTAL COST. S 
VESSEL MATERIAL COST. $ 

TOTAL MATERIALS COST FOR TES. J 
TOTAL FAPR1CATI ON COST FOR TES* * 
TES TOTAL COST. $ 

VFSSEL TOTAL COST. 1 
.* FABRICATION COSTS ARE SCALED FRO* 
TUBE MATERIAL COST COEFFICIENT 
FELIX DIAMETER. M. *.F10.6) 


MR )TE !6» 
MR ITE 46. 

106 ) 
107 ) 

TSHR 

MR 1TE (6. 

1C6) 

L SHE LL 

MR 1TE (6. 

111) 

FSHE LL 

MR 1TE (6. 

109 ) 

to 

MR 1TE (6. 

11C ) 

01 

MR ITE ! 6 . 

I 12 ) 

XLT 

MR ITE 46. 

1 13 ) 

XNT 

MP ITE 1 6 • 

791 ) 

OHEL 

MR 1 T E ( 6 . 

114) 

TINS 

MR ITE (6. 

136) 

000 

MR ITE (6 . 

139 ) 

0 1 0 

MR ITE (6. 

1 15 ) 

7 SHE LL 

MP )TE (6. 

116) 

TEND 

MR I TE ( 6 . 

I 17) 

RHOM 

MR ITE (6. 

1 IB ) 

RHOST 

MP ITE 4 6 • 

1 19) 

RHOF 

MRITE46. 

132 ) 

RH01 NS 

MR ITE (6. 

120 ) 

T MC INS 

MR JTE(6. 

12 1 ) 

THCST 

MR 1TE(6 . 

122 ) 

THCF 

MR ITE (6. 

123) 

FRF 

MR ITE 46 . 

124) 

V1SCF 

MR ITE 46. 

125 ) 

TAMP 

MR ITE 46. 

129 ) 

TSHP 

MR ITE 46. 

126 ) 

UNIT 

MR ITE 46. 

127) 

CPST 

MR ITE 46 . 

135 ) 

► TVE SS 

MR ITE 4 6. 

136) 

MTTUEE 

MP ITE 46 . 

137) 

VTINS 

MR ITE46. 

133 ) 

FGTT43T 

MR ITE 4 6 » 

134 ) 

MTTOT 

MR ITE 46. 

770 ) 

ASPECT 


MCf 


• « F 10.4) 
*.F 10.4 ) 
*.F10.4) 


'.FI 2.2) 
*.F 12.2) 
*.F 12.2 ) 
* . F 12.2) 
*.F 12.2) 
*.F 12.2) 
*.F 12.2) 

* . F 1 2 .2 I 
•.F12.2) 
*.F 12.2) 
*.F 12.2) 
• . F 12.?) 
*.F 12.2) 
*.F 1? .2) 
••F 12.2) 

• * F 1 2 .2 ) 
• * F 1 £ .2 ) 

I TES*./) 

*.F 1 C .4 ) 
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I F ( J5 . EC. . 1 ) GC TO 4 
Wfi 1 TE ( 6 • 1 26 ) XL T 

4 *F ITEC6.750 ) 

WF ITE ( 6 .75 1 ) FMATS 
WFITE(6,752) FFSSA 
KFITE(6.7E3) FFSSW 
WF ITE (6 .754 ) FFSTA 
WFITE(6.755) FFSTW 
WRITEt6.7E6) FFSEA 
WF 1 TE ( 6 .757 ) PF5BV 
tt F I TE ( 6 .758 ) GATOPS 
WF ITE (6.759 } CASI DF 
WRITE (6.760) CABOTS 
WF ITE (6.761 ) FMyW 
WFITE(6,762) FMVE 
WF I TE (6.763 ) FFVV 
W F ITE (6.764 ) FFVW 
tt F 1 TE ( 6 • 790 ) CPMTW 
WF1TE(6.765) FMTW 
WF ITE (6.766 ) FFTV 
KFITE(6.767) I F T A 
tt F ITE (6.760 ) F TK 
WF ITE (6.769) F INS 
ttF ITE (6 .769 ) 
tt F ITE (6.771 ) 

W F 1 Tt « 6 * 77? ) CMTUPE 
W F 1 TE (6.773) COTUPE 
WF ITE (6.774 ) CCSTOR 
WF1TE(6. 775 ) COINS 
tt F ITE (6.776 ) CNSIDE 
W F I TE ( 6 • 777 ) COSIDE 
tt F ITE(6»778) CMBOTS 
HP ITE (6,779) COBOTS 
WR ITE (6.760 ) CMTCPS 
HP ITE (6*761 ) COTOPS 
tt F 1 TE (6.762) CMSHRD 
WF ITE(6.7B3) COSH FD 
ttFlTE(6.764) CMVCSS 
WRITE(6.7eb) CCVESS 
ttF ITE(6 .765) CMTOT 
tf P ITE ( 6 .786 ) CCFTES 
WRITE(6»767) COTOT 

t INITIALIZE TAELE OF TLEPENT ENERGIES AT T 1 N I T 

CALL T X PROP (TINIT .ENEL.FRAC.1.1 .J6) 

E- 0. 

DC 5 1 = 1. NE 

E T ( 1 )=T INI T 
2 ( 1 )=FR AC 
E =E + ENER ( I ) 

5 CONTINUE 

ttF ITE ( 6 . 194 ) E 
K 1 —2 

C DECIDE FLOW OF NO -FLEW 

2 IF ( K 7. NE • 1 ) GC TO 15 
OF L = 0 
GC TO 20 

C CALCULATE INSIDE HEAT TRANSFER COEFFICIENT AND RESISTANCE 

15 RE=WF*CFE 
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450 1 F (RE .LT .RE CR1T) GO TC 451 
46? HT C1=CHTC*RES*0.P 

GO TO 463 

451 HT Cl=3.fc6*THCF/DI 

RE M0D=RE4 (Dl/rHEL )="*C.5 
QX=(<REM0D**2)4PRF )$$0.25 

ZNU=0.181«CX*1 1 .-0.e30/GX + 35.4/GX*<=2-207./OX<*34 4 19/QX**4 ) 
HTCI=2NU4HTC 1 
453 RF L 1= 1 ./t HTC I$AT UBE 1 ) 

C COMPUTATION OF PUPP MQU1REMENT FOR PUMP EFF=1.C 

PD = • 1 02SFF ( RE )#XLT*(XF/FLAREA)*<2./<2.*D1*RHCF ) 

IF CPD.GT.PDMAX) PDMAX=PD 
PCWER=WF*PD/C 1 02 • 4RHCF ) 

IF (POWER. GT .P) AX ) PMAX^POMER 
*DRK=P0WER*D1 1ME/36CC. 

#PUMP = WPUMP-»WC'R< 

IF (WF.GT.WPEAK ) WPEAK = *F 
2 C CONTINUE 

66 IFIK7.EQ.2) Gt TO 13C 
C ENERGY CALCULATIONS - COOLING OR NO FLOW 

TF ( 1 )=TL 
DO 50 1=1. NE 
IF (K7.EQ . I 1 GO TO 2 1 
IF C Z f 1 l.EQ.O.) R S OL L = 0 
IF tZ( 1 1 .EQ.O. ) GC TC 66 
V S T<=Z C I )*ELM./(RH0M*XNT 1 
DSTK=SQRT < 1 .273*VSTK$XNE/XLT4DC «*2 ) 

RSOLL = ALOGC D ST K/D 01/ ITHCSt T ( 1 ) ) <CRK TL ) 

66 RTDT = RFL URTUEEL + RSOLL 
XN TU=1 ./CPTOT=C=WF#CPF ) 

IF CXNTU.GT. IOC.) CFF=1. 

IF (XNTU.GT . IOC.) GO TO 140 
EF F = 1 .-EXPI-XNTU) 

140 TF(I*1)=TFII)-»EFF*(T< 1 )-TF ( I ) ) 

GF L=WF*CPF*(TF <1 )-TF I 1 + 1 ) ) 

6 DECIDE TOP < 7 1 ) • B0TTCM(72). OR MIDDLL 

21 CONTINUE 
QT OT=QFL 

ENERI 1 )=ENERC l)4QT0T<n 1ME 

CALL TK PROP (TCI) . ENERfl ) .2 ( I ) • 1 .2 • J6 ) 

GO TD 50 

< BOTTOM ELEMENT CALCULATIONS 

71 RA XLW = 1 . /< 2 ./RTUPEA-* 2./RSHELA+2 ./RNA ( 1 1 )4RA>E 

RA XUP=1 .✓ ( 1 ./RTUPEA4 l./RSHELA+l ./ ( (RMAC 1 M RM A C 1 4 1 ) ) 22 . ) ) 

CAXLW = UAME-TCI1 l/RAXLK 

GAXUP=(Tf 1-»1)-T*1 ) ) / RAX UP 

QT0T=GAXLW4QA>UP+GFL*CINS 

ENER( I ) = ENER II )4QT0T*CT I ME 

CALL TKPROPCTUl.ENEPU ).2 (I ).2 J6) 

GAXLW=-GAXUP 
60 TO 50 

C TIP ELEMENT CALCULATIONS 

72 RA XUP=1 ./(2./RTUBEA4 2./RSHELA42 ./RMA(I 1 )4RAXE 
QA XUP=t TAME-T I 1 ) l/PAXUP 
QT0T=GAXLW4QA>UP4G1N£-»0FL 

ENER{ l) = ENER(IMOTOT=fDTlME 

CALL TK PR OP C T ( I ) * ENE R 1 1 1 .2 ( I 1 .2 .2 . J6 1 

*0 continue 
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GL TO 193 

C ENERGY CALCULAT JONS - HEATING 

13C TF ( NE 4 1 ) = TU 

DO 150 1=1. NE 
Jr NE — 1 + 1 

IF (Z< Jl.EQ.O.) RSCLL=C 
IF CZC JJ.EO.O.) GO TO 176 
VSTK=Z(J)*EL*/(RHOM*>NT) 

DSTK=SGRT ( 1 .273* VST X *XNE /X L T 4D0 <*2 ) 

R£ OLL=ALOG( DSTK/D O) / ( T HCS ( T ( J) > *CRX TL) 

17C RT0T=RFLI+RTUEEL4RSCLL 
XNTU=1 ./(RTOT**F*CPF } 

IF (XNTU.GT.1CC.) EF F = 1 » 

IF CXNTU.GT. 100.) GO TO 1AI 
EF F=1 .-f XP(-XNTU) 

1 A 1 TF ( J)=TF( J41 HEFF«(T (J)-TF < J4 1 ) > 

GF L = WF*CPF*(TF F J-» 1 > - TF ( J » ) 

QT OT=GFL 

ENER( J)=ENER ( J)*CTOT =50 T I ME 

CALL TXPPOPIT(J) .FNEM J).Z ( Jl.l .2.J6) 

GC TO lfO 

( UP ELEMENT CALCULAT ICNS 

171 R A XUP= 1 ./(2./RTU6EA42./RSHELA+2 ./RMACJ) ) 4 F AXE 

RAXLW = 1 ./< 1 ./RTUPEA-t 1 ./RSHELA4 1 ./ ( (RVA( JHGPA { J- 1 ) ) /2 . ) ) 

GAXUP=(TAWE-T( J) )/RAXUP 

GAXLW=(T( J-l l-TC J ))/RAXLX 

QT0T=GAXUP4CAXL*401NS4QF L 

ENERC J)=ENER < J ) 4 0 TO T SOT I ME 

CALL TK PROP C T ( J) »ENER(J)»Z (J).2.2.J6) 

C A XUP = ~GA XL Vi 
GC TO 150 

C BOTTOM ELEMENT CALCULATIONS 

172 R A XLW= 1 ./( 2./RTUPF A + 2./RSHELA+2 . /RM A tJ))-*RA>E 
GAXLW=(TAMB“T (J) l/RAXLM 

Q T OT = QAXL *4Q AXUP4 G IN S4QF L 

EKER( J)=ENER I j>4GTCT=tniME 

CALL T K PROP ( T t J) ,ENEF(J).Z(J).2.2*J6) 

ICC CONTINUE 

193 IF(CT.E0.T5) CO TO 1 ( t 
IF (CT.EG.T6I GO TO ICC 
CT 1=CT1 4DT IME 

IF (CT.GE.T5.AFD.CT.LE.Tfc) C FGT M = 1 C T -T5 ) /3fc 00 . 

IF (CT1.LT.160C.) GO TO 257 
CCMT=(WF$TF (NE 4 1 >4»fCCLL«TL )/(WF4WC0LL> 
fc‘RITE(6» 24 0 I f T . C HGT F • W'F , T F (NE 4 1 ) .WCOLL .CCMT 
CT 1=0. 

257 RETURN 

165 K R 1 TE (6.165) 15 
CT 1=0. 

GO TO 190 

166 *RITE(6.166) It 
CPGTM = ( T6-T5 1/36 00. 

ISO E = 0 . 

W R 1 TE (6.187) 

IF (K7.EC. 1 J GO TO 168 
DT 167 1=1, NL 
T E MPEL= 1 • 8* T (1 1-459.67 

X R 1 TE ( 6 . 188 ) 1 .T( I ) ,Z( 1 ) .TFC I ) . TEMPEL 
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1 C 7 E = E+ENER ( 1 ) 

L C WE ND = NE ■* 1 

W R 1 TE (6.1 92 ) LCWE fJD . 7f CNE-» 1 ) 

GC TD ieo 

168 D C. 169 1=1* NE 

Wfi HE (6*189) ]»T(11*Z(I) 

169 E=E+ENER(I) 

160 KP 17E (6.191 ) L 

WR17E(6*24C) E 7 • C HOT »' . WF • 7 F (NE ■* 1 ) • WCOLL • CGM7 
CHGTM=0. 

198 RETURN 
END 
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SUBROUTINE TNTROP (TE MP.ENTH .FRA C .M.K.2. J6 ) 

C THIS SUBROUTINE ESTABLISHES TABLES 01 E l t RENT ENT1ALPY VS. THE 

r FRACTION OF SOLID THERMKEEP AND TEMPERATURE 

C 0 MMON/PROP 7 K/ELM .CPST.STHM.STME 
DIMENSION EM(14).EE(14),T!14).H(14).Z!14) 

2CC F CRMAT ! 12 ) 

ZC1 FORMAT! F10.2.F12.S.F 10.6) 

202 FORMAT! 1-H1.5X.* TABLE OF TEMPERATURE. ENTHALPY. AND SOLID FRACTION 
* OF THERMKEEF • ./✓) 

2C3 F GRMAT ! 1 OX . * T!I)*.1CX,» H!I)».10X.» 2!J)*./) 

2 04 F C RMA T (tX.F 10.2.6X.F 12.5.5X.F10.6) 

20b FORMAT ( 5X • • ENERGY BEYOND TABLE LIMITS FOR ENT H= * .E 10 .4 ) 

GO TO! 1 .2 ) *tC2 

1 IF ! J6.EG.2) GO TO 1C 1 
MR ITEI6.202 ) 

RF AD! 5.200 I K4 
MR I TE 16*203) 

101 J4=K4-1 
J3=l 

DO 11 1 = 1. IC4 

IF ! J6.EC.2 ) GO TO 102 

R E AD 15.201 ) TlD.Ht 1 J.Zll) 

1C2 CONTINUE 

IF 1K2.EG.2 ) GO TO 2 

IF ! J3.EC.2) 60 TO 10 

IF IT ( I l.LT -TEMP) GO TO 10 

ENTH=ELM*! H< 1-1 >4 (TEMP-T ! I - 1 ) ) * !H ( I ) -H ! I - 1 ) )✓ ! T C 1 )-T ! 1 -1 ) ) ) 

FRAC=Z! D + lTEMP-TlI-l ))*!Z! I )-Z ! 1 - 1 ) )/! T 1 1 >-T ( I - l ) ) 

J3 =2 

10 EM(I|=ELM*H( 1 )4STMM4CPST«T ! I ) 

EE ! 1)=ELM*H! I )4ST ME 4CFS74T ( 1 ) 

IF ! J6.EQ.2) GO TO 1 1 
MRI7E16.2G4) T ! 1 ) .H ( 1 ) • Z ! I > 

11 CONTINUE 
RE TURN 

2 GO TO! 3.4) .M 

3 IF (ENTH.GT.EM1K4 ) .OR .ENTH.LT.EM 11 )) GO TO 2&r 
DC 5 1=1. J4 

IF IENTH.GT. EM! 1 + 1 ) ) CO TO 5 

T E MP=T 1 I H IT! Ml ) -T U ) )#!E NTH-E M( I ) )/!EN< 1+] ) - £ A <I) ) 

FR AC = 2 ! I )4 12 * I-»l | -2 f I ) >*!E NTH-E.M * 1 » )/«EM I l-*l )-EN « 1 > ) 

RE TURN 
f CONTINUE 
GO TO 2E>0 

4 I F 1ENTH.GT.EE !K4 ) .OR ,f NTH. LT. EE (1 ) ) GO TO 2b C 
DC 6 1=1. J4 

IF (ENTH.GT .EE ! M 1 )) GO TO 'j 

T F MP=T 1 I)4!T1J41)-T!1))*!E NTH-EE ( I ) ) / ! EE 1 1 4 1 ) -E E ( 1 ) ) 

F R AC = Z ! I )4 !Z ! 1 4 1 ) -2 ! I ))*( E NTH-EE ! I) )/!Et 1 14 1 )-f E 1 1 ) ) 

RE TURN 
t CONTINUE 

2&C MRITE16.205) F NTH 
RE TURN 
END 
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FUNCTION THCS(T) 

C THIS FUNCTION COMPUTFS THE THERMAL CONDUCTIVITY CF SOLID THFRMKEEP 

C AT 250F. IT IS NOW A DUMMY FUNCTION WITH A CONSTANT VALUE 

C IN KJ/M-SEC-N 

tf cs=i .oeooE-: j 

RETURN 

END 
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FUNCTION THCLO) 

C THIS FUNCTION COMPUTES THE THERMAL CONDUCTIVITY LF LIQUID THPHCEEP 

C AT 750F. IT 15. NC»» A DUMMY FUNCTION WITH A CLNSlAKT VALUE 

C IN KJ/M-SEC-K 

T F CL= 1 • 0207E -3 
RE TURN 
END 
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FINCTION F F ( R ) 

TFIS FUNCTION COKFUTLS THE CIRCULAR TUBE FRICTION FACTOR F OR 
ROUGHNESS RATIO OF .CC012 
IF (R.GT .2000. ) GO TO 1 
FF =64 ./R 
RE TURN 

1 IF (R.GT. 500000. I GO TO 2 
FF =.33e7*F5*«(-.2E49 ) 

RE TURN 
? F F =0 .0124 
RE TURN 
END 
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